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We study the merger of black hole- neutron star binaries with a variety of black hole spins aligned or 
antialigned with the orbital angular momentum, and with the mass ratio in the range Mbh/Mns = 
2-5, where Mbh and Mns are the mass of the black hole and neutron star, respectively. We model 
neutron-star matter by systematically parametrized piecewise polytropic equations of state. The 
initial condition is computed in the puncture framework adopting an isolated horizon framework 
to estimate the black hole spin and assuming an irrotational velocity field for the fluid inside the 
neutron star. Dynamical simulations are performed in full general relativity by an adaptive-mesh 
refinement code, SACRA. The treatment of hydrodynamic equations and estimation of the disk mass 
are improved. We find that the neutron star is tidally disrupted irrespective of the mass ratio when 
the black hole has a moderately large prograde spin, whereas only binaries with low mass ratios, 
Mbh/Mns < 3, or small compactnesses of the neutron stars bring the tidal disruption when the 
black hole spin is zero or retrograde. The mass of the remnant disk is accordingly large as > O.IMq, 
which is required by central engines of short gamma-ray bursts, if the black hole spin is prograde. 
Information of the tidal disruption is reflected in a clear relation between the compactness of the 
neutron star and an appropriately defined "cutoff frequency" in the gravitational-wave spectrum, 
above which the spectrum damps exponentially. We find that the tidal disruption of the neutron star 
and excitation of the quasinormal mode of the remnant black hole occur in a compatible manner in 
high mass-ratio binaries with the prograde black hole spin. The correlation between the compactness 
and the cutoff frequency still holds for such cases. It is also suggested by extrapolation that the 
merger of an extremely spinning black hole and an irrotational neutron star binary does not lead to 
the formation of an overspinning black hole. 

PACS numbers: 04.25. D-, 04.30.-w, 04.40.Dg 



I. INTRODUCTION 

Coalescing binaries composed of a black hole (BH) 
and/or a neutron star (NS) are among the most promis- 
ing sources of gravitational waves for ground-based 
laser-interferometric gravitational-wave detectors such 
as LIGO [J and VIRGO 0. Detections of gravita- 
tional waves will be accomplished in a decade to come 
by planned next-generation detectors such as advanced 
LIGO, advanced VIRGO, and LCGT Q. Because gravi- 
tational waves are much more transparent to the absorp- 
tion and scattering by the material than electromagnetic 
waves and even neutrinos are, gravitational- wave astron- 
omy is expected to become a powerful and unique way 
to observe strongly gravitating phenomena in our Uni- 
verse. Among such phenomena, the merger of a BH-NS 
binary plays an important role to investigate properties 
of the NS such as the radius and the equation of state 
(EOS) of a high-density nuclear matter [J-Q ■ An impor- 
tant constraint on the EOS is obtained from detection of 
a 1.97 ± 0.04M^ NS, which is the most massive NS cur- 
rently known [10| . by recent pulsar-timing observation. 
However, we still do not know the realistic EOS of the 
NS because there is no robust measurement of the NS 
radius 11]. To constrain the NS radius and EOS by ob- 



serving gravitational waves from the BH-NS binary, we 



have to prepare accurate theoretical templates of gravi- 
tational waveforms employing a wide variety of the NS 
EOSs and other physical parameters. For this purpose, 
numerical relativity is the unique approach. 

The merger of a BH-NS binary is also an important 
target for the astrophysical study because it is a poten- 
tial candidate for the progenitor of short-hard gamma- 
ray bursts (GRBs) in the so-called merger scenario (see 
Refs. [TU, [l3| and references therein for reviews) . If a NS 
is tidally disrupted during the merger of a BH-NS binary, 
a system composed of a spinning BH and a hot, massive 
accretion disk of > O.OIMq may be formed. Such an 
outcome could be a central engine of the GRB, because 
it could radiate a large amount of energy > 10 48 erg via 
neutrino emission or the Blandford-Znajek mechanism 
[bij in a short time scale of < 2 s and hence could launch 
a GRB jet. The merger scenario of the GRB is attractive 
when a short-hard GRB is associated with a galaxy of low 
star-formation rate [HI, [l6[ because the collapsar model 
of GRBs [13] is not preferable. Only numerical relativ- 
ity can answer quantitatively the question whether the 
formation of a massive accretion disk is possible in the 
BH-NS binary merger. 

Fully general relativistic study of BH-NS binaries has 
achieved progress in recent years, both in computations 
of quasiequilibrium states [18l422l | and in dynamical sim- 
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ulations of the merger [23l - [3^ |. Although these include 
several simulations of spinning BH-NS binaries with a 
qualitative T = 2 ideal-gas EOS, to date, only limited 
number of simulations have been performed taking into 
account both the nuclear-theory based EOS and the BH 
spin [30l ]. In particular, we still do not understand the 
dependence of the merger process and resulting gravita- 
tional waveforms on the BH spin and the EOS of the NS 
in detail. One goal of current numerical relativity is to 
clarify the effect of the BH spin on the BH-NS binary 
merger and subsequent NS tidal disruption adopting a 
wide variety of the NS EOSs. 

In this paper, we report our latest results obtained by 
numerical-relativity simulations with a variety of the NS 
EOSs and the BH spins. We employ five piecewise poly- 
tropic EOSs (see Sec. IIIBj) . all of which do not conflict 
with the current observation of the 1.97 ± 0.04Af Q NS. 
We systematically choose physical parameters such as BH 
mass, NS mass, and BH spin in an astrophysically real- 
istic range. While we onlyconsider relatively low mass 
BHs in the previous work [3l| , we adopt a wider range of 
the BH mass because a large BH spin enhances the NS 
tidal disruption for high mass-ratio binaries. We clarify 
the dependence on the BH spin and on the NS EOS of 
the properties of the merger remnants and characteristics 
of gravitational waves — in particular, the gravitational- 
wave spectrum. 

This paper is organized as follows. In Sec. HH we de- 
scribe methods for a solution of initial conditions, piece- 
wise polytropic EOSs, and models of BH-NS binaries em- 
ployed in this paper. In Sec. IIII1 the formulation and 
methods of numerical simulations are summarized. Sec- 
tion IIVI presents the numerical results and clarifies the 
effect of the BH spin and NS EOS on the tidal disrup- 
tion, merger remnants, and gravitational waves. Section 
fVl is devoted to a summary. Throughout this paper, we 
adopt the geometrical units in which G = c = 1, where 
G and c are the gravitational constant and the speed 
of light, respectively. Our convention of notation for 
physically important quantities is summarized in Table 
|TJ The nondimensional spin parameter of the BH, to- 
tal mass of the system at infinite separation, mass ratio, 
and compactness of the NS are defined as a = Sbh/-^bh> 
m = M BH + M NS , Q = Mbh/M ns , and C = M NS /R NS , 
respectively. Latin and Greek indices denote spatial and 
spacetime components, respectively. 



II. INITIAL CONDITION 

As in our previous works [1^, [3l[ , we employ BH-NS 
binaries in quasiequilibrium states for initial conditions 
of our numerical simulations. In this section, we summa- 
rize the formulation and methods for the computation 
of a quasiequilibrium state, specifically our method of 
estimating the spin angular momentum of the BH in a 
binary and of determining the position of the rotation 
axis. The details of the formulation and numerical meth- 



TABLE I. Our convention of notation for physically important 



quantities. 


Symbol 


M irr 


mi ■ l "i l r 1 1 n t t 

the irreducible mass 01 the BH 


•5BH 


Ine magnitude 01 tne BH spin angular momentum 


Mbh 


The gravitational mass of the BH in isolation 


Mns 


The gravitational mass of the NS in isolation 


-Rns 


The circumferential radius of the NS in isolation 


Mo 


The Arnowitt-Deser-Misner mass of the system 


m 


The total mass of the system at infinite separation 


Q 


The mass ratio 


c 


The compactness of the NS 


a 


The nondimensional spin parameter of the BH 



ods, except for the issues on the BH spin, are described 
in Ref. [221 ] . to which the reader may refer. Computa- 
tions of the quasiequilibrium states are performed using 
the spectral- method library LORENE |33| . 

A. Formulation and methods 

We compute a quasiequilibrium state of the BH-NS bi- 
nary as a solution of the initial value problem of general 
relativity As far as the orbital separation is large 

enough, the time scale of the orbital contraction due to 
the gravitational radiation reaction is much longer than 
the orbital period, and, therefore, we may safely neglect 
the gravitational radiation reaction in the calculation of 
the quasiequilibrium state. In numerical simulations of 
the binary coalescences, we have to track > 5 orbits in 
order to calculate accurate gravitational waveforms dur- 
ing the late inspiral and merger phases, and hence, the 
orbital separation of the initial condition has to be large 
enough. For such initial conditions, we can neglect the 
gravitational radiation reaction. Thus, we give a BH- 
NS binary in a quasicircular orbit, i.e., the binary in an 
approximate equilibrium state in the corotating frame. 
To satisfy the quasiequilibrium requirements described 
above, we assume the existence of a helical Killing vector 
with the orbital angular velocity f2, 

e = (d t r + n(d v y. (i) 

We also assume that the NS is in the hydrostatic equilib- 
rium in the corotating frame and has an irrotational ve- 
locity field, which is believed to be a reliable approxima- 
tion to an astrophysically realistic configuration [35l [36j . 

We compute the three- metric 7^ , the extrinsic curva- 
ture Kij, the lapse function a, and the shift vector p l 
by a mixture of the extended conformal thin-sandwich 
approach [37], [38l and the conformal transverse-traceless 
decomposition|34| in the puncture framework (3*j - |4l| . In 
this formalism, we assume the conformal flatness of the 
three- metric 7^ = ip 4 7ij — ^fiji the stationarity of the 
conformal three- metric £^ij = 0, the maximal slicing 
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condition K = r ) lJ Kij = 0, and its preservation in time 
£^K = 0, where /y is the flat spatial metric. Assuming 
that the puncture is located at x P , we set the conformal 
factor ip and a weighted lapse function $ = atp as 



1> 



1 



M P 
2r B H 



+ <?!>, $ = 1 



2r B H 



where Mp and M$ are positive constants of mass dimen- 
sion, and ?'bh = \x % — Xp\ is a coordinate distance from 
the puncture. We adjust Mp to obtain a desired mass 
of the BH, Mbh, and determine M$ so as to satisfy the 
virial relation, i.e., the equality of the Arnowitt-Deser- 
Misner (ADM) mass and the Komar mass, which holds 
when the spacetime is stationary and asymptotically flat 
[42l [43l ] . (j), f3 l , and 77 are determined by solving ellip- 
tic equations derived from the Hamiltonian constraint, 
momentum constraint, and quasiequilibrium conditions, 
£^ij = and £^K = 0. We note that these quasiequi- 
librium conditions can be replaced by 9t7y = and 
dt K = in a conformal flatness approximation. In the 
puncture framework, we also decompose a conformally 
weighted extrinsic curvature = ip 2 Kij as 



K 



-p 



(3) 



where Wi is an auxiliary three- vector field and V* is the 
covariant derivative associated with . Kfj is a singular 
part of the extrinsic curvature, which is associated with 
the linear and spin angular momenta of the BH [iij . 



i.P 



BH 



l 3 P 



BH 



k nBHi 



{fij hlj)l Pk 



BH 



— 3 — [ e kilSpl lj + €)-jlSpl 



(4) 



BH 



where l l = z BH /rBH is a unit radial vector, k = fijl 3 , 
and eijk is the Levi-Civita tensor associated with the flat 
metric fij. Pf n and Sp are parameters associated with 
the linear and spin angular momenta of the BH, respec- 
tively. Here we determine -P BH by the condition in which 
the total linear momentum of the binary vanishes, and 
we adjust Sp to obtain a desired spin angular momentum 
of the BH, S'gjj. The elliptic equation to determine Wi 
is obtained by taking the derivative of Eq. ([3]) and using 
the momentum constraint. 

The spin angular momentum of the BH, 5 B h: is evalu- 
ated on the apparent horizon (hereafter AH), S, accord- 
ing to the isolated horizon framework (see Refs. [45|, |46| 
for reviews) . Because we do not know the position of the 
AH in advance in the puncture framework, we have to 
determine the location of the AH numerically [47[- On 
the numerically determined AH, an approximate rota- 
tional Killing vector may be defined using the method 
developed by Cook and Whiting [48( with the normaliza- 
tion condition proposed by Lovelace and his collaborators 
[49l | . We focus only on the case in which the BH spin is 
aligned or antialigncd with the angular momentum of the 
binary in this work, and hence, the axis of the BH spin 



is uniquely determined. Hereafter we set the rotational 
axis of the binary, equivalently the axis of the BH spin, to 
be the z axis and consider only the approximate Killing 
vector (j) 1 associated with the rotation in this direction. 
Using <p l , we obtain the spin angular momentum of the 



(2) BH Sgjj = S'gjj via the surface integral at the AH, 



5bh " 8^ 



1 



(5) 



We adjust Sp to obtain a desired value of 5"bh (hereafter 
Sbh)- We note that Sp and S^ do not agree exactly in 
the BH-NS binary spacetime due to the contribution to 
the extrinsic curvature from the NS, associated with Wi- 
Because we adopt a conformal flatness approximation 
for the induced metric, the Christodoulou mass of the 
BH evaluated on the AH, M H = y/M^ + (^ H /4M^ r ), 
and the gravitational mass evaluated at spatial infinity, 
Mbh; do not agree even for a single BH system due to the 
presence of so-called junk waves. This difference leads to 
an ambiguity in defining the nondimcnsional spin param- 
eter of the BH. Here, we define the nondimensional spin 
parameter of the BH with respect to the mass evaluated 
at spatial infinity, i.e., 



5bh 



(0) 



The reason for this is that the mass and nondimensional 
spin parameter of the BH evaluated at the AH quickly 
(in our simulations, within ~ 1 ms) relax to Mbh and 
a, defined at spatial infinity, respectively, as the BH ab- 
sorbs the junk radiation in the vicinity of the BH (49l[50j. 
We note that these values show the damping oscillation 
before the relaxation in the same manner as the "scalar- 
curvature spin" of Ref. [49| shows, because our method of 
evaluating these values in the simulation is basically the 
same as the method to define the scalar-curvature spin 
in Ref. ji| (see Scc.|HTB|). 

To compute the equations of hydrostatic equilibrium 
for the NS matter, we assume an ideal fluid, for which 
the energy-momentum tensor is described by 



T^ v = phu^u v + Ptf 



(7) 



where p is the rest-mass density, P is the pressure, h = 
l + e+ (P/p) is the specific enthalpy, e is the specific 
internal energy, and is the four-velocity of the fluid. 
Basic equations for the hydrostatics are derived from the 
condition of irrotation, or the vanishing of the vorticity 
two-form dH-il, 



V M (/iu I/ ) - Vvihup) = 0, 



(8) 



and the conservation of the specific momentum of the 
fluid along the helical Killing vector field £^(hu lJ ) = 0. 
The specific enthalpy is determined from the first integral 
of the rclativistic Eulcr equation (rclativistic Bernoulli in- 
tegral), and other thcrmodynamical quantities are sub- 
sequently obtained by using the EOS, which is described 
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in Sec. Ill Bl The four-velocity of the fluid is calculated 
using the velocity potential "J, as hu,i = in the as- 
sumption of the irrotational velocity field. The elliptic 
equation for VP is derived from the equation of continu- 
ity, V^pu") = 0. 

We have no definite condition to determine the loca- 
tion of the center of mass of the binary in the puncture 
framework and we use this ambiguity to reduce an un- 
physical initial orbital eccentricity. We found in our pre- 
vious work (22], [2t| that orbits with a small eccentricity 
could be obtained using the "3PN-J method," i.e., a phe- 
nomenological method to determine the location of the 
rotational axis in which the total angular momentum of 
the binary for a given value of fimo agrees with that 
calculated from the third post-Newtonian (3PN) approx- 
imation. In the present case in which the BH has a finite 
amount of spin angular momentum, we extend the previ- 
ous 3PN-J method to include the contribution from the 
BH spin up to the 2.5PN order [55|, [56|. Specifically, the 
location of the rotational axis is chosen from the condi- 
tion that the orbital angular momentum of the binary 
agrees with a sum of nonspin terms given by Eq. (4 ) of 
Ref. [57| and of spin terms given by Eq. (7.10) of Ref. |56j 
for a given value of Omo- We show in Sec. IIV Al that this 
extended 3PN-J method again leads to the initial condi- 
tion with a small eccentricity. 



B. Piecewise polytropic equations of state 

The matter inside the NS in the late inspiral phase is 
believed to be well-approximated by a zero-temperature 
nuclear matter because the cooling time scale of the NS 
in typical BH-NS binaries is shorter than the time scale of 
the gravitational radiation reaction (58j . Hence, we em- 
ploy a cold EOS, for which the rest-mass density, p, deter- 
mines all other thermodynamical quantities, for calculat- 
ing the quasiequilibrium state of the BH-NS binary. To 
model nuclear-theory-based EOSs at high density with 
a small number of parameters, we employ a piecewise 
polytropic EOS. It is a phcnomcnologically parametrized 
EOS of the form 



P(p) = Kip Tl for pi_i < p < pi (l<i< n), 



(9) 



where n is the number of the pieces used to parametrize 
an EOS, pi is the rest-mass density at the boundary 
of two neighboring ith and (i + l)th pieces, is the 
polytropic constant for the ith piece, and Ps is the adi- 
abatic index for the ith piece. Here, po — 0, p n — > 00, 
and other parameters (pi,Ki,Ti) are freely chosen. Re- 
quiring the continuity of the pressure at each pi, 2n 
free parameters — say (Ki,Ti) — determine the EOS com- 
pletely. The specific internal energy, e, and hence the spe- 
cific enthalpy, h, are determined by the first law of ther- 
modynamics and continuity of each variable at boundary 
densities, pi. 

It was shown that piecewise polytropic EOSs with four 
pieces approximately reproduce most properties of the 



nuclear-theory-bascd EOSs at high density |59j ■ If we fo- 
cus on low mass NSs with relatively low central density, 
the EOS at high density plays a minor role. Thus, we 
adopt a simplified piecewise polytropic EOS composed 
of two pieces, one of which models the crust EOS and 
the other of which the core EOS. This simplification is 
based on the fact that NSs in the observed binary NSs 
often have fairly small masses < 1.4M Q [6(| and the max- 
imum rest-mass density in such NSs may not be so high 
that the EOS at high density plays only a minor role in 
determining their structure. Furthermore, the maximum 
rest-mass density inside the NS should only decrease dur- 
ing the evolution of the BH-NS binary due to tidal elon- 
gation of the NS by the companion BH. 

Tabic HIl lists the EOSs which we employ in this study. 
Following Refs. 0, [3l|, we always fix the EOS for the 
crust region by the parameters below: 



Ti = 1.35692395, 
ki/c 2 = 3.99873692 x 1(T 8 (g/cm 3 ) 1 ^ 1 



(10) 
(11) 



The EOS of the core region is determined by two param- 
eters. One is the adiabatic index of the core EOS, P2- 
Whereas we find that properties of gravitational waves 
and the merger remnants depend on T2 in our previ- 
ous study [3l| . we always fix T2 — 3 in this work to 
focus on clarifying the effect of the BH spin aside from 
the difference in the adiabatic index. The other param- 
eter is chosen to be the pressure p at a fiducial density 
Pfidu = 10 14,7 g/cm 3 because p is closely related to the 
radius and deformability of the NS [61J]. We vary the 
value of p systematically to investigate the effect of the 
stiffness of the EOS [62| . With the given values of T2 and 
p, K2 and p\ are determined as 



K-2 



Pi 



PPfilu 



(KlA2) 1/(IWl) . 



(12) 
(13) 



C. Models 

Numerical simulations are performed for a wide range 
of nondimcnsional BH spin parameter, a, as well as for a 
variety of the mass ratios, Q. For nonspinning BH-NS bi- 
naries, we already found that the low mass ratio of Q < 3 
is required for tidal disruption of NSs to occur sufficiently 
outside the innermost stable circular orbit (ISCO) of the 
BH unless the EOS is extremely stiff [H|]. If the tidal 
disruption occurs inside or at an orbit very close to the 
ISCO, we do not see strong effects of the tidal disrup- 
tion. In such cases, gravitational waveforms are similar 
to those of a BH-BH binary even in the merger phase, and 
the mass of the remnant disk is negligible [63| . However, 
the allowed range of the mass ratio for the tidal disrup- 
tion is modified drastically for a BH-NS binary with the 
prograde BH spin f64l . |65| because the ISCO radius [66| 
of the BH with a prograde spin becomes smaller by a 
factor of 1 - 6 67J than that of the nonspinning BH with 
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TABLE II. Key ingredients of the adopted EOSs. T2(~ 3.0) is the adiabatic index in the core region and p is the pressure at the 
fiducial density psdu = 10 14 ' 7 g/cm , which determines the polytropic constant K2 of the core region and pi, the critical rest-mass 
density separating the crust and core regions. M max is the maximum mass of the spherical NS for a given EOS. -Ri3s(-Ri2, ^145) 
and Ci35(Ci2, C145) are the circumferential radius and the compactness of the NS with Mns = 1.35Mq(1.2M0, I.45M0). 



Model 


r 2 ioi 


r ,ioP (g/cm 3 ) 


pi(10 14 g/cm 3 ) 


M max [M ] 


-R135 (km) 


C135 


R12 (km) 


Cu 


Rias (km) 


C145 


2H 


3.0 


13.95 


0.7033 


2.835 


15.23 


0.1309 


15.12 


0.1172 


15.28 


0.1401 


1.5H 


3.0 


13.75 


0.9308 


2.525 


13.69 


0.1456 


13.63 


0.1300 


13.72 


0.1561 


H 


3.0 


13.55 


1.232 


2.249 


12.27 


0.1624 


12.25 


0.1447 


12.27 


0.1744 


HB 


3.0 


13.45 


1.417 


2.122 


11.61 


0.1718 


11.60 


0.1527 


11.59 


0.1848 


B 


3.0 


13.35 


1.630 


2.003 


10.96 


0.1819 


10.98 


0.1614 


10.93 


0.1960 



the same mass. Strong spin effects for the tidal disrup- 
tion are also found in the numerical-relativity simulation 
of the spinning BH-NS binary merger with a simplified, 
T-law EOS [28l ] . In this paper we perform a more sys- 
tematic study of the tidal disruption for different EOSs, 
masses of each component, and BH spins. 

Table Hill summarizes several key quantities for the ini- 
tial conditions in our numerical simulations. The label 
for the model denotes the EOS name, the mass ratio, 
the NS mass, and the nondimensional spin parameter 
of the BH. Specifically, "a75," "a5," and "a-5" corre- 
spond to the spin parameters a = 0.75, 0.5, and —0.5, 
respectively. For example, HB-Q3M135a5 means that 
the EOS is HB and (Q,M NS ,a) = (3, 1.35M Q , 0.5). Al- 
though we vary the NS mass systematically, the results 
of the merger remnant are reported only for binaries with 
Mns = 1.35 A/© in this paper because the difference in 
the NS mass complicates the properties of the remnant, 
such as the mass of the disk. Results for Mns 7^ 1.35M 
are analyzed only for gravitational waves. 

For the same value of the mass ratio, we basically pre- 
pare the initial conditions with the same value of the 
initial angular velocity f2o normalized by the total mass 
of the binary, f2o m o- For 2H EOS, in which the NS radius 
is the largest, we exceptionally adopt a smaller value of 
Qomo than for other EOSs to guarantee > 5 orbits before 
tidal disruption occurs. The reason for this is that the 
tidal disruption occurs for a large orbital separation in 
2H EOS. When the BH has a prograde spin, the number 
of orbits to the merger for a given value of Homo increases 
due to spin-orbit repulsive interaction [68j], compared to 
the nonspinning BH case. On the other hand, when the 
BH has a retrograde spin, the number of orbits decreases 
due to spin-orbit attractive interaction. For a = —0.5, 
the number of orbits is typically by ~ 1 orbit smaller 
than for a = 0. For this reason, we also prepare the ini- 
tial condition with a smaller value of ^o m o for H EOS 
and a = —0.5. 

III. METHODS OF SIMULATIONS 

Numerical simulations are performed using an 
adaptive-mesh refinement (AMR) code SACRA [&§] . The 
formulation, the gauge conditions, the numerical scheme, 



and the methods of dia gnos tics are basically the same as 
those described in Ref. [31(, except for the correction in 
the treatment of hydrodynamic equations in a far region. 
Thus, we here only briefly review them and describe the 
present setup of the computational domain for the AMR 
algorithm and grid resolution. 

A. Formulation and numerical methods 

SACRA solves the Einstein evolution equations in the 
Baumgarte-Shapiro-Shibata-Nakamura (BSSN) formal- 
ism [20, [7lJ with the moving- puncture gauge f39j - l4il ] . It 
evolves a conformal factor W = 7 -1 ^ 6 , the conformal 
metric 7^ = 7 _1 ' 3 7ij, the trace of the extrinsic curva- 
ture K, a conformally weighted trace- free part of the ex- 
trinsic curvature Aij = 7" 1/,3 (-?Qj — Kjij), and an aux- 
iliary variable T l = —dj^ 13 . Introducing an auxiliary 
variable B l and a parameter ?/ s , which we typically set 
to be ~ Mbr/Mq in units of c = G = M © = 1, we 
employ a moving-puncture gauge in the form [72j 

(St- P j dj)a = -2aK, (14) 
{d t -^d j )p i = ^/A)B\ (15) 

(d t - pdAB* = (d t - pdj)r - Vs b\ (i6) 

We evaluate the spatial derivative by a fourth-order cen- 
tral finite difference, except for the advection terms, 
which arc evaluated by a fourth-order noncentered, up- 
wind finite difference, and employ a fourth-order Rungc- 
Kutta method for the time evolution. 

To solve the hydrodynamic equations, we evolve p* = 
pau t W~ 3 , Ui = hui, and = hau 1 — P / (pau 1 ) . The ad- 
vection terms are handled with a high-resolution central 
scheme by Kurganov and Tadmor [73| with a third-order 
piecewise parabolic interpolation for the cell reconstruc- 
tion. For the EOS, we decompose the pressure and spe- 
cific internal energy into cold and thermal parts as 

P = Pco\d + Pth , £ = fcold + £th- (17) 

We calculate the cold parts of both variables using the 
piecewise polytropic EOS from the primitive variable p, 
and then the thermal part of the specific internal energy 
is defined from e as e t h = £ ~ £ C o\d- Because £th vanishes 
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TABLE III. Key parameters and quantities for the initial conditions adopted in numerical simulations. The adopted EOS, mass 
ratio (Q), NS mass in isolation (Mns), nondimensional spin parameter of the BH (a), initial angular velocity (Qo) in units of 
c 3 /Gmo, baryon rest mass (M»), compactness of the NS in isolation (C), maximum rest-mass density (p ma x), ADM mass of the 
system (Mo), and total angular momentum of the system (Jo), respectively. See also pTO ] for models of nonspinning BH-NS 
binaries. 



Model 


EOS 


Q 


M ns [Mq\ 


a 


Gfl m /c 3 


M*[M Q 


C 


Pmax (10 g/Cm ) 


71 J r 71 JT 1 

M [Mq\ 


J [GMq/c\ 


2H-Q2M135a75 


2H 


2 


1.35 


0.75 


0.025 


1.455 


0.1309 


3.740 


4.014 


13.83 


1.5H-Q2M135a75 


1.5H 


2 


1.35 


0.75 


0.028 


1.468 


0.1456 


5.104 


4.012 


13.42 


H-Q2M135a75 


H 


2 


1.35 


0.75 


0.028 


1.484 


0.1624 


7.019 


4.012 


13.42 


HB-Q2M135a75 


HB 


2 


1.35 


0.75 


0.028 


1.493 


0.1718 


8.263 


4.012 


13.42 


B-Q2M135a75 


B 


2 


1.35 


0.75 


0.028 


1.503 


0.1819 


9.762 


4.012 


13.42 


2H-Q2M135a5 


2H 


2 


1.35 


0.5 


0.025 


1.455 


0.1309 


3.740 


4.014 


14.02 


1.5H-Q2M135a5 


1.5H 


2 


1.35 


0.5 


0.028 


1.468 


0.1456 


5.104 


4.012 


13.63 


H-Q2M135a5 


H 


2 


1.35 


0.5 


0.028 


1.484 


0.1624 


7.018 


4.012 


13.63 


HB-Q2M135a5 


HB 


2 


1.35 


0.5 


0.028 


1.493 


0.1718 


8.263 


4.012 


13.63 


B-Q2M135a5 


B 


2 


1.35 


0.5 


0.028 


1.503 


0.1819 


9.762 


4.012 


13.63 


2H-Q2M135a-5 


2H 


2 


1.35 


-0.5 


0.022 


1.455 


0.1309 


3.740 


4.019 


15.15 


H-Q2M135a-5 


H 


2 


1.35 


-0.5 


0.025 


1.484 


0.1624 


7.018 


4.017 


14.74 


HB-Q2M135a-5 


HB 


2 


1.35 


-0.5 


0.028 


1.493 


0.1718 


8.262 


4.015 


14.41 


B-Q2M135a-5 


B 


2 


1.35 


-0.5 


0.028 


1.503 


0.1819 


9.760 


4.015 


14.41 


2H-Q2M12a75 


2H 


2 


1.2 


0.75 


0.025 


1.282 


0.1172 


3.465 


3.568 


10.93 


H-Q2M12a75 


H 


2 


1.2 


0.75 


0.028 


1.303 


0.1447 


6.421 


3.566 


10.60 


HB-Q2M12a75 


HB 


2 


1.2 


0.75 


0.028 


1.310 


0.1527 


7.523 


3.566 


10.60 


B-Q2M12a75 


B 


2 


1.2 


0.75 


0.028 


1.317 


0.1614 


8.833 


3.566 


10.60 


2H-Q2M145a75 


2H 


2 


1.45 


0.75 


0.025 


1.572 


0.1401 


3.926 


4.312 


15.96 


H-Q2M145a75 


H 


2 


1.45 


0.75 


0.028 


1.607 


0.1744 


7.452 


4.309 


15.48 


HB-Q2M145a75 


HB 


2 


1.45 


0.75 


0.028 


1.617 


0.1848 


8.811 


4.309 


15.48 


B-Q2M145a75 


B 


2 


1.45 


0.75 


0.028 


1.629 


0.1960 


10.46 


4.309 


15.48 


2H-Q3M135a75 


2H 


3 


1.35 


0.75 


0.028 


1.455 


0.1309 


3.737 


5.357 


20.00 


1.5H-Q3M135a75 


1.5H 


3 


1.35 


0.75 


0.030 


1.468 


0.1456 


5.100 


5.355 


19.64 


H-Q3M135a75 


H 


3 


1.35 


0.75 


0.030 


1.484 


0.1624 


7.013 


5.355 


19.64 


HB-Q3M135a75 


HB 


3 


1.35 


0.75 


0.030 


1.493 


0.1718 


8.256 


5.355 


19.64 


B-Q3M135a75 


B 


3 


1.35 


0.75 


0.030 


1.503 


0.1819 


9.753 


5.355 


19.63 


2H-Q3M135a5 


2H 


3 


1.35 


0.5 


0.028 


1.455 


0.1309 


3.737 


5.357 


20.36 


1.5H-Q3M135a5 


1.5H 


3 


1.35 


0.5 


0.030 


1.468 


0.1456 


5.100 


5.356 


20.02 


H-Q3M135a5 


H 


3 


1.35 


0.5 


0.030 


1.484 


0.1624 


7.012 


5.356 


20.01 


HB-Q3M135a5 


HB 


3 


1.35 


0.5 


0.030 


1.493 


0.1718 


8.255 


5.356 


20.01 


B-Q3M135a5 


B 


3 


1.35 


0.5 


0.030 


1.503 


0.1819 


9.753 


5.356 


20.01 


HB-Q3M135a-5 


HB 


3 


1.35 


-0.5 


0.030 


1.493 


0.1718 


8.253 


5.359 


21.46 


2H-Q3M145a75 


2H 


3 


1.45 


0.75 


0.028 


1.572 


0.1401 


3.923 


5.754 


23.07 


H-Q3M145a75 


H 


3 


1.45 


0.75 


0.030 


1.607 


0.1744 


7.445 


5.751 


22.65 


HB-Q3M145a75 


HB 


3 


1.45 


0.75 


0.030 


1.617 


0.1848 


8.803 


5.751 


22.65 


B-Q3M145a75 


B 


3 


1.45 


0.75 


0.030 


1.629 


0.1960 


10.45 


5.751 


22.65 


2H-Q4M135a75 


2H 


4 


1.35 


0.75 


0.030 


1.455 


0.1309 


3.735 


6.702 


26.07 


H-Q4M135a75 


H 


4 


1.35 


0.75 


0.032 


1.484 


0.1624 


7.007 


6.700 


25.62 


HB-Q4M135a75 


HB 


4 


1.35 


0.75 


0.032 


1.493 


0.1718 


8.249 


6.700 


25.63 


B-Q4M135a75 


B 


4 


1.35 


0.75 


0.032 


1.503 


0.1819 


9.746 


6.700 


25.62 


2H-Q4M135a5 


2H 


4 


1.35 


0.5 


0.035 


1.455 


0.1309 


3.732 


6.698 


25.64 


H-Q4M135a5 


H 


4 


1.35 


0.5 


0.035 


1.484 


0.1624 


7.004 


6.698 


25.63 


HB-Q4M135a5 


HB 


4 


1.35 


0.5 


0.035 


1.493 


0.1718 


8.244 


6.698 


25.63 


B-Q4M135a5 


B 


4 


1.35 


0.5 


0.035 


1.503 


0.1819 


9.740 


6.698 


25.63 


2H-Q5M135a75 


2H 


5 


1.35 


0.75 


0.036 


1.455 


0.1309 


3.730 


8.044 


30.95 


H-Q5M135a75 


H 


5 


1.35 


0.75 


0.036 


1.484 


0.1624 


7.000 


8.044 


30.95 


HB-Q5M135a75 


HB 


5 


1.35 


0.75 


0.036 


1.493 


0.1718 


8.241 


8.044 


30.95 


B-Q5M135a75 


B 


5 


1.35 


0.75 


0.036 


1.503 


0.1819 


9.736 


8.043 


30.95 
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in the absence of shock heating, eth is regarded as the 
finite-temperature part. In this paper, we adopt a T-law 
ideal-gas EOS for the thermal part, 

Pth = (r th - l)pe th , (18) 

to determine the thermal part of the pressure, and choose 
r t h equal to the adiabatic index in the crust region, Fi, 
for simplicity. 

Because the vacuum is not allowed in any conservative 
hydrodynamic scheme, we put an artificial atmosphere 
of a small density outside the NS in the same way as 
done in our previous work [3l|. The total rest mass of 
the atmosphere is always less than 10 _4 Mq, and hence, 
we can safely neglect spurious effects by accretion of the 
atmosphere onto the remnant disk as long as the disk 
mass is much larger than 10 _ -M©. 



B. Diagnostics 

We extract gravitational waves by calculating the out- 
going part of the Weyl scalar \&4 at finite coordinate radii 
r = 400-800A/© and by integrating ^4 twice in time as 

h+(t) - ih x (t) = - f dt' J dt"V A (t"). (19) 

In our previous works [29l l3l| , we directly perform this 
integration of ^(t) and then subtract a quadratic func- 
tion of the form ait 2 + a\t + ao to eliminate unphysi- 
cal drift components in the waveform, using the least- 
square fitting to obtain constants ao, 01, and a 2 . In this 
work, we adopt a "fixed-frequency integration" method 



proposed by Reisswig and Pollney [74| to obtain gravi- 
tational waveforms with less unphysical components. In 
this method, we first perform a Fourier transformation 
of ^4 as 

* 4 H = J dm A (t)e iut . (20) 
Using this, Eq. p^|) is rewritten as 

h+(t) - ih x (t) = — — f ^le-^du. (21) 

Z7T J UJ £ 

We then replace 1/lu 2 of the integrand with 1/ujq for 
I a; I < ujq, where ujo is a positive free parameter in 
this method. By appropriately choosing luq, this proce- 
dure suppresses unphysical, low-frequency components of 
gravitational waves. As proposed in Ref. [74j |. we choose 
loq to be ~ O.SmQ for m 7^ mode gravitational waves, 
where m is the azimuthal quantum number. For the 
m = mode gravitational waves, we adopt u>o ~ O.8f2o 
and confirm that our results depend only very weakly 
on this choice. We also adopt this method to calculate 
the energy AE and angular momentum A J radiated by 
gravitational waves. Exceptionally, we adopt the previ- 
ous method of direct time integration to estimate the or- 
bital eccentricity in the inspiral phase because the fixed- 
frequency integration method may change the modula- 
tion in the gravitational waveform. 

For comparisons between numerically calculated grav- 
itational waveforms and those calculated in the PN ap- 
proximations, we use the Taylor-T4 formula for two-point 
masses in circular orbits |75| with an additional contri- 
bution from the BH spin angular momentum (76j . In this 
formula, the time evolution of the orbital angular veloc- 
ity £l(t) and orbital phase Q(t) are computed using a 
nondimensional angular velocity X(t) = [m Sl(i)] 2 / 3 by 



dX MvX z 



dt 5mo 



, / 743 11 \ „ / 113 19 , v3/2 

1 1 v ) X + 4tt y H va X i/2 

336 4 J I 12 A 6 



/34103 r 2 13661 59 2 \ v2 

h 5y 2 -I v H v X 2 

V 18144 A 2016 18 / 

/4159 189 \ /31571 1165 \ 3 , /21863 79 ,\ 1 . /2 

! v 7r + v Y + -X ~ v v a > X b > 2 

\ 672 8 / V 1008 24 j A 4 A V !008 6 / J 

f 16447322263 1712 16 , / 56198689 451 9 \ 541 , 5605 , 

-< 7F H 7T 7T V -\ V 1/ 

\ 139708800 105 )iu 3 \ 217728 48 j 896 2592 
856, ,„„__,. 80tt /64153 457 \ 2 /20 1135 \ 1 - 

/4415 358675 91495 ,\ /2529407 845827 41551 ,\ ,„ 2 

V V 2 7T + V H V 2 \ - 127TY 2 

\4032 6048 1512 / \ 27216 6048 864 / A A 

/1505 v\ - / 1580239 451597 2 2045 , 107 a \ 1 v7/ J 

f — ) y — ( v 2 H v z H vx a \X 7 ' 2 , (22 

V 24 8 J A V 54432 6048 432 6 A J J J v ' 

I 

dO _ X 3 / 2 where v = Q/(l + Q) 2 , X = aQ/{l + Q), and 7jE « 0.5772 

dt mo 
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is the Euler constant. After X(t) and 0(i) are obtained, 
we calculate the complex gravitational-wave amplitude 



h 22 of the (l,m) = (2,2) mode and the spectrum up to 
the 3PN order using the formula shown in Refs. (76l.l77||. 
Here, h 22 is 



/7TW10 2i0 

" 8 V5~ e X 
2173 1069 



107 
l2~ 



2047 



1512 216 1512 
f 27027409 856 2 2 
1 646800 ~ 105 7£; + 3^ 



X 2 - 



107 34 



V 7T 

21 21 



428 
105 



ln(16A) 



24wj A 5/2 
( 278185 41 ,\ 20261 

7T V 1 

V 33264 96 / 2772 



+ ■ 



114635 
99792 1 



428 
105' 



(24) 



where D is the distance between the center of mass of 
the binary and an observer. Hereafter, we simply refer 
to this formula as the Taylor-T4 formula irrespective of 
the presence of the BH spin. Another way for deriving 
an approximate waveform is to employ an effective one- 
body approach (see Ref. [78[ and references therein for 
reviews). In accompanied papers (63|, comparisons be- 
tween numerical waveforms and those of the effective one 
body approach are extensively performed. 

To estimate the mass of the remnant disk, we calculate 
the total rest mass outside the AH 

M r>rAH = I p*d 3 x, (25) 

Jr>r A H 

where tah = !'ah(^ 'p) is the radius of the AH as a func- 
tion of the angular coordinates. We note that we system- 
atically underestimated disk masses in our previous works 
performed with an old version of SACRA [2^, [3l| , because 
we evolved hydrodynamic variables and estimated disk 
masses only in the finer domains (described in Sec. IIII C[) 
of the size ~ 200 3 km 3 . Such a domain size is insufficient 
for the estimation of the disk mass if tidal disruption oc- 
curs at a distant orbit, especially for the case in which 
the NS radius is large (~ 15 km). In this study, we cor- 
rect the treatment of hydrodynamics and the estimation 
of disk masses: We follow the hydrodynamics for a wide 
computational domain of the size 1000 3 -2000 3 km 3 . We 
still possibly underestimated disk masses because some 
of the material escapes from our computational domains 
and we cannot follow their return which would occur if 
they are bounded. 

We determine key quantities of the remnant BH, i.e., 
the mass Mbhj and nondimcnsional spin parameter at, 
from the circumferential radius of the AH, assuming that 
the deviation from the Kerr spacetime is negligible in the 
vicinity of a BH horizon. We estimate the remnant BH 
mass, MBH,f , from the circumferential radius of the AH 
along the equatorial plane C e divided by An, i.e., C e /47r, 
which gives the BH mass in the stationary vacuum BH 
spacetime. Similarly, the nondimcnsional spin parameter 



of the remnant BH, at, is estimated from the ratio of the 
circumferential radius of the AH along the meridional 
plane C p to C e using the relation 

This also holds for the stationary vacuum BH with the 
nondimensional spin parameter af. Here, f + = 1 + 
\Jl — a 2 is a normalized radius of the horizon, and E(z) 
is an elliptic integral 

r /2 , 

E(z)= / VI - zsin 2 0d6. (27) 

For comparison, the nondimcnsional spin parameter of 
the remnant BH is also estimated from C e and the irre- 
ducible mass of the remnant BH Mj rri f using the relation 

M irr , f = -^Jl + ^l-a^, (28) 

which holds for the stationary vacuum BH. The spin pa- 
rameter obtained using this relation is referred to as a{2 
according to Ref. [2^|. Finally, we also estimate af from 
the values of the remnant BH computed using approxi- 
mate conservation laws 

M BH ,c = M - M r>rAH - AE, (29) 

^BH,c = ^0 ~ Jr>r A H ~ (30) 

where the total angular momentum of the material lo- 
cated outside the AH, J r >r AH , is approximately defined 
by 

Jr>rAH — / p-fhUtpd? 'X. (31) 
Jr>r A H 

Here, we assume that the orbital angular momentum of 
the BH is negligible. The nondimcnsional spin parame- 
ter of the remnant BH is defined by = Jbb.c/Mbh c 
again according to Ref. [29j . 
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C. Setup of AMR grids 



In SACRA, an AMR algorithm is implemented so that 
both the radii of compact objects in the near zone and 
the characteristic gravitational wavelengths in the wave 
zone can be covered with sufficient grid resolutions simul- 
taneously. Our AMR grids consist of a number of compu- 
tational domains, each of which has the uniform, vertex- 
centered Cartesian grids with (2N + 1, 2N + 1, N + 1) 
grid points for (x, y, z) with the equatorial plane symme- 
try at z = 0. We always choose N = 50 for the best 
resolved runs in this work. We also perform simulations 
with N = 36 and 42 for several arbitrary chosen models 
to check the convergence of the results and find approx- 
imately the same level of convergence as that found in 
the previous work (see the Appendix of Ref. [3lj]). I n 
the Appendix „ of this paper, we show the convergence 
of gravitational waveforms and the masses of the rem- 
nant disks. The AMR grids are classified into two cat- 
egories: one is a coarser domain, which covers a wide 
region, including both the BH and NS, with its origin 
fixed at the approximate center of mass throughout the 
simulation. The other is a finer domain, two sets of 
which comove with compact objects and cover the re- 
gion in the vicinity of these objects. We denote the edge 
length of the largest domain, the number of the coarser 
domains, and the number of the finer domains by 2L, l c , 
and 21 f, respectively. Namely, the total number of the 
domains is l c + 2lf. The grid spacing for each domain is 
hi = L/(2 l N), where I = 0-(l c + If - 1) is the depth of 
each domain. 

Table IIVI summarizes the parameters of the grid struc- 
ture for our simulations. The structure of the AMR grids 
depends primarily on the mass ratio of the binary because 
the distances between two objects and the center of mass 
depend strongly on the mass ratio for our initial mod- 
els. Specifically, we choose (l c , If) — (4, 4) for all binaries 
with M NS = 1.35M and Q = 2,3, and 4. We choose 
(l c , If) = (3, 5) for binaries with Q = 5. For binaries with 
^-/ns 7^ 1.35M©, we choose (l c ,lf) = (3, 4) because we do 
not evaluate disk masses for them. In all the simulations, 
L is chosen to be larger than or comparable to the grav- 
itational wavelengths at an initial instant Ao = 7r/^o- 
One of the two finest regions covers the semimajor axis 
of the NS by ~ 42-45 grid points. The other covers the 
coordinate radius of the AH typically by ~ 20 grid points, 
depending on the mass ratio and the BH spin. For the 
Q = 5 runs, the total memory required is about 11 G 
bytes. We perform numerical simulations with personal 
computers of 12 G bytes memory and of core-i7X proces- 
sors with clock speeds of 3.2 or 3.33 GHz. We use 2-6 
processors to perform one job with an OPEN-MP library. 
The typical computational time required to perform one 
simulation (for ~ 50 ms in physical time of coalescence 
for the a = 0.75 case) is 4 weeks for the 6 processor case. 



IV. NUMERICAL RESULTS 

We present numerical results of our simulations, focus- 
ing, in particular, on their dependence on the BH spin 
and NS EOS. First, we review general merger processes 
in Sec. [TO] Sections HVBl HVCl and|lVD]are devoted 
to the analysis of properties of the remnant disk and BH 
formed after the merger. Gravitational waveforms are 
shown in Sec. IIVE1 their spectra in Sec. IIVF1 and the 
energy and angular momentum radiated by gravitational 
waves in Sec. IIVG1 



A. Overview of the merger process 

Figure [T] plots the evolution of the coordinate sepa- 
ration defined by x l sep = x l NS — x l BU for models HB- 
Q2M135a5, HB-Q2M135, and HB-Q2M135a-5, for which 
J]omo takes the same values. Here, x^ s is the position of 
the maximum rest-mass density and x Bil is the location 
of the puncture, x p . Figure Q] shows that the numbers of 
orbits increases as the BH spin increases from retrograde 
to prograde [28| . Specifically, the number of orbit are 
~ 7, 5.5, and 4 for a = 0.5, 0, and —0.5, respectively. 
This difference comes primarily from the spin-orbit in- 
teraction between these two angular momenta ff^ : in 
the PN approximation, a force proportional to the in- 
ner product of the orbital and spin angular momenta of 
two objects appears at 1.5PN order. Here, we do not 
have to consider the NS spin angular momentum in the 
assumption of the irrotational velocity field and, there- 
fore, we only consider the interaction between the orbital 
and BH spin angular momenta throughout this paper. 
When these two angular momenta are parallel and the 
inner product is positive (a > 0), an additional repul- 
sive force works between the BH and NS. This repulsive 
force reduces the orbital angular velocity because the cen- 
trifugal force associated with the orbital motion can be 
reduced, and hence, the luminosity of gravitational ra- 
diation, which is proportional to f2 10 ' 3 , is also reduced. 
This strong dependence of the luminosity on O makes the 
approaching velocity smaller in the late inspiral phase, 
and, therefore, the number of orbits increases. Con- 
versely, when these two angular momenta are antiparallcl 
(a < 0), an additional attractive force increases the an- 
gular velocity and gravitational-wave luminosity in the 
late inspiral phase. In this case, the orbital separation 
decreases faster due to a larger approaching velocity, and 
the number of orbits becomes smaller as the retrograde 
BH spin increases. All these results agree qualitatively 
with those of Ref. [Hj]. 

The fate of BH-NS binaries is classified into two cate- 
gories. One is the case in which the NS is disrupted by 
the BH tidal field before the BH swallows the NS, and the 
other is the case in which the BH swallows the NS with- 
out tidal disruption. In this paper, we focus mainly on 
the former case. We plot snapshots of the rest-mass den- 
sity profiles and the location of the AH on the equatorial 
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TABLE IV. Setup of the grid structure for the simulation with our AMR algorithm. l c and If are the number of coarser domains 
and a half of finer domains, respectively. Ax = hi = L/(2 l N) (I = l c + lf — 1) is the grid spacing at the finest-resolution domain 
with L being the location of the outer boundaries along each axis. iidiam/Ax denotes the grid number assigned inside the 
semimajor diameter of the NS. Ao is the gravitational wavelength of the initial configuration. (See Ref. [3l[ for models with 
nonspinning BHs.) 



Model 


lc 


h 


Ax/Mo 


R diam /Ax 


L/Xo 


2H-Q2M135a75 


4 


4 


0.0471 


90.8 


2.386 


1.5H-Q2M135a75 


4 


4 


0.0426 


87.7 


2.417 


H-Q2M135a75 


4 


4 


0.0377 


86.2 


2.138 


HB-Q2M135a75 


4 


4 


0.0347 


87.1 


1.968 


B-Q2M135a75 


4 


4 


0.0324 


86.7 


1.837 


2H-Q2M135a5 


4 


4 


0.0471 


90.8 


2.378 


1.5H-Q2M135a5 


4 


4 


0.0426 


87.7 


2.410 


H-Q2M135a5 


4 


4 


0.0377 


86.2 


2.131 


HB-Q2M135a5 


4 


4 


0.0347 


87.2 


1.962 


B-Q2M135a5 


4 


4 


0.0324 


86.7 


1.831 


2H-Q2M135a-5 


4 


4 


0.0470 


90.7 


2.092 


H-Q2M135a-5 


4 


4 


0.0376 


86.4 


1.902 


HB-Q2M135a-5 


4 


4 


0.0347 


87.1 


1.962 


B-Q2M135a-5 


4 


4 


0.0324 


86.7 


1.831 


2H-Q2M12a75 


3 


4 


0.0583 


84.7 


1.476 


H-Q2M12a75 


3 


4 


0.0442 


85.3 


1.252 


HB-Q2M12a75 


3 


4 


0.0410 


85.7 


1.162 


B-Q2M12a75 


3 


4 


0.0389 


84.2 


1.102 


2H-Q2M145a75 


3 


4 


0.0461 


85.2 


1.166 


H-Q2M145a75 


3 


4 


0.0347 


85.3 


0.985 


HB-Q2M145a75 
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FIG. 2. Evolution of the rest-mass density profile in units of g/cm 3 and the location of the AH on the equatorial plane for 
model HB-Q3M135a75. The filled circle denotes the region inside the AH. The color (gradational) panel on the right of each 
plot show log 10 (p). 
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FIG. 3. The same as Fig. H but for model HB-Q3M135a5. 
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FIG. 4. The same as Fig. [2] but for model HB-Q3M135a-5. 



plane at selected time slices for models HB-Q3M135a75, 
HB-Q3M135a5, and HB-Q3M135a-5 in Figs.UHl respec- 
tively. The NS is disrupted outside the IS CO in the a > 
cases (Figs. [2] and [3]) and forms a one-armed spiral struc- 
ture with a large angular momentum. The material in 
the inner part of the spiral arm gradually falls onto the 
BH due to angular momentum transport via hydrody- 
namic torque in the spiral arm. The material with a 
sufficiently large specific angular momentum escapes the 
capture by the BH and forms an accretion disk, which 
survives for a time much longer than the dynamical time 



scale ~ a few ms. We note that the prompt infall of the 
one-armed spiral structure onto the BH occurs from a rel- 
atively narrower region for a = 0.5 than for a = 0.75. The 
reason is that the inner edge of the spiral arm contacts 
the AH well before the arm becomes nearly axisymmetric 
due to a large radius of the AH and IS CO for a = 0.5. 
The infall of the disrupted material from a narrow region 
of the BH frequently occurs when the NS is tidally dis- 
rupted in a binary with a high mass ratio, whereas this 
is rare in a binary with a nonspinning BH because the 
NS is not disrupted in a high mass-ratio binary. This 
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difference in the merger process is well-reflected in grav- 
itational waveforms (see Sec. llVE]) . By contrast, the NS 
is swallowed by the BH without tidal disruption, and es- 
sentially no material is left outside the ISCO for model 
HB-Q3M135a-5 (Fig.g]). 

Note that the feature of the NS tidal disruption ap- 
pears very weakly not only for model HB-Q3M135a-5 
but also for model HB-Q3M135 (a = 0) because the mass 
ratio Q = 3 is so high that the tidal effect is less impor- 
tant for the nonspinning BH with the typical NS radius 
~ 11-12 km. The enhancement of the tidal effect by a 
prograde BH spin results primarily from the decrease of 
the BH ISCO radius In the Boyer-Lindquist coordi- 
nates, a Kerr BH has an ISCO with a smaller radius than 
a Schwarzschild BH by a factor of > 1/6, depending on 
a for a prograde orbit: The ISCO radius approximately 
halves when the BH spin increases from a — to 0.75. 
On the other hand, the orbital separation at the onset of 
mass shedding depends only weakly on the BH spin in 
the Boyer-Lindquist coordinates |80h82| |. This decrease 
of the ISCO radius enhances the possibility for the dis- 
rupted material to escape capture by the BH and to form 
a more massive remnant disk than in the nonspinning BH 
case. The retrograde BH spin plays an opposite role; the 
ISCO radius of the Kerr BH increases by a factor of 1—1.5 
for a retrograde orbit, and hence, the tidal effect is less 
important in the merger process. 

Before closing this subsection, we estimate the degree 
of (undesired) orbital eccentricity in our simulations to 
assess the circularity of the orbital motion. For this pur- 
pose, we compute the evolution of the gauge-invariant 
orbital angular velocity Q(t), which is defined from the 
(I, m) = (2, 2) mode of #4 by 



fi(f) 



1 



2 | / $ 4 (l =m = 2)dt\ ' 



(32) 



The evolution of the orbital angular velocity in our sim- 
ulation agrees with that derived from the Taylor- T4 for- 
mula in the inspiral phase within a small modulation, 
typically AQ/f2 < 5%, which is equivalent to the orbital 
eccentricity of < 3%. This amount of orbital eccentricity 
is as small as that observed in the nonspinning BH case 
with a low mass ratio Q = 2 [3l| . 



B. Global properties of the disk 

The mass of the remnant disk reflects the significance 
of NS tidal disruption in a clear way because the disk 
formation is a result of tidal disruption. A massive disk 
is formed if tidal disruption of the NS occurs far outside 
the ISCO. If the mass shedding starts in the vicinity of or 
inside the ISCO, only a small portion of the mass is left 
outside the AH. The material is not left outside the AH 
when the mass shedding does not occur before the BH 
swallows the NS, and the merger of a BH-NS binary may 
be indistinguishable from that of a BH-BH binary except 
for very small tidal corrections to the inspiral. Thus, the 



mass of a remnant disk is a reliable indicator of the degree 
of tidal disruption. 

Figure [5] plots the time evolution of the rest mass lo- 
cated outside the AH, M r>rAli , for Q = 2 and 3 with 
different nondimcnsional BH spin parameters a — 0.75, 
0.5, 0, and -0.5. In both plots, M NS = 1.35M Q and 
HB EOS are adopted. We note that the results revised 
from the previous work [3~i| are plotted for a = 0. The 
dependence of M r>rAll on a for HB EOS found here is 
similar to those for other EOSs. We set the time origin 
to be an approximate merger time t merger . These plots 
indicate that the mass of the material left outside the AH 



relaxes to a quasisteady value for t — t T 



> 3-4 ms, 



and the relaxed value increases monotonically as the BH 
spin increases from retrograde to prograde. This is con- 
sistent with the decrease of the BH ISCO radius with the 
increase of its spin, as described in Sec. IIV Al In partic- 
ular, the remnant disk mass at « 10 ms after the merger 
is > 0.1M Q for all the EOSs with (Q,a) = (2, > 0.5) 
and (< 4, 0.75), as shown in Tabled and > 0.05M Q for 
(Q,a) = (3,0.5), irrespective of the EOS. The formation 
of such a massive disk may be encouraging for the BH- 
NS binary merger hypothesis of a short-hard GRB. For 
the a = —0.5 cases, by contrast, massive accretion disks 
of > 0.01M© are not expected to be formed as merger 
remnants even for Q = 2 unless the EOS is extremely 
stiff (the NS radius is ~ 15 km). This fact indicates that 
the retrograde BH spin is unfavorable for producing a 
central engine of a short-hard GRB. 

The prograde BH spin enhances the disk formation 
dramatically for a BH-NS binary with a high mass ra- 
tio, for which the disk mass is very low when the BH is 
nonspinning. We plot the time evolution of M r>rAll for 
Q = 4 and 5 with different EOSs in Fig. [SI In both plots, 
M~ns = l-35Af Q and a = 0.75 are adopted. Figure [5] 
clearly shows that a massive accretion disk is formed for 
Q = 4 and 5 if the BH has a prograde spin of a = 0.75. 
Namely, the formation of a massive accretion disk is uni- 
versal for the merger of a BH-NS binary with a mass 
ratio of Q < 5 as far as a ~ 0.75 and Mns = 1.35M 
(cquivalcntly, Mbh ^ 6.75M Q ). Note that a heavy BH 
of Mbh ^ 5M Q is predicted to be realistic as an astro- 
physical consequence of the stellar evolution with solar 
metallicity [H| (see, e.g., Ref. [H| for a population syn- 
thesis study) and hence as a possible progenitor of the 
short-hard GRB. 

For more quantitative discussion, we plot the disk mass 
estimated at « 10 ms after the merger for all the models 
with Q = 2 and for models with (Q, a) = (3, > 0) as a 
function of the NS compactness, C, in Fig. [71 Numerical 
values of M r>rAH are shown in Table fVl as well as other 
quantities associated with the merger remnants. For any 
fixed value of a, a negative correlation between M r>rAll 
and C is found to hold in Fig. [71 This correlation indicates 
that the NS with a larger compactness is less subject 
to tidal deformation and disruption than the NS with 
a smaller compactness for any fixed value of a. This 
correlation is expected from the nature of a tidal force as 
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FIG. 5. Evolution of the rest mass of the material located outside an AH of the BH, M r>rAH . The left and right panels show 
the results for Q = 2 and 3, respectively. In both plots, Mns = 1.35ATq and HB EOS are adopted. The results revised from 
the previous paper [3lJ are plotted for a = 0. 




FIG. 6. The same as Fig. [5]but for different models Q = 4 (left) and 5 (right). In both plots, Mns = 1.35Mq and a = 0.75 are 
adopted. 
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a finite-size effect, as found in the study of nonspinning 



BH-NS binaries [31|. On the other hand, Fig. [7] again 
shows that the progradc BH spin increases the disk mass 
for any fixed value of C. A remarkable fact is that the disk 
mass does not decrease steeply to a value of <C 0.1M Q 
as the compactness increases for binaries with (Q,a) = 
(< 3, > 0.5). We expect that the coalescence of a BH-NS 
binary with (Q, a) = (< 3, > 0.5) may always produce 
a remnant disk of > 0.01M Q within a plausible range 
of the NS compactness, C < 0.2, although it is possible 
only if C < 0.18 for (Q,a) = (2,0) and C < 0.16 for 
(Q,a) = (2,-0.5) or (3,0). 

The dependence of the disk mass on the NS compact- 
ness is different for different values of the mass ratio. 
We plot in Fig. [8] the disk mass as a function of the 
NS compactness as in Fig. but for a = 0.75 and 
0.5. This figure shows that the disk mass depends more 
strongly on C when the mass ratio, Q, is larger. The 
disk mass is larger for smaller values of Q when the EOS 
is soft and C > 0.16, except for HB-Q2M135a75 and 
HB-Q3M135a75, for which the disk masses depend only 
weakly on Q. This dependence on Q is expected from 
the comparison between the mass-shedding radius, r s hod, 
and the ISCO radius, nsco, 



nsco 



-2/3 



(33) 



where we assume Newtonian gravity for simplicity. This 
relation states that a larger amount of mass can escape 
the capture by the BH and can form an accretion disk 
when Q is small because the mass shedding sets in at 
relatively more distant orbit. However, the disk mass 
may be larger for larger values of Q when the EOS is 
stiff as C < 0.15 for a > 0.5 and 2 < Q < 5. This 
should be ascribed to the redistribution process of the 
specific angular momentum of the NS to the disrupted 
material and to subsequent behavior of the material (such 
as collision of the fluid elements in spiral arms). This 
feature suggests that a binary with a larger value of Q, 
sa Y Q si 6, possibly form a massive remnant disk of > 
O.IMq and could be a progenitor of a short- hard GRB if 
the EOS is stiff and the BH has a large spin > 0.5. 

To clarify the dependence of the disk mass on the BH 
spin, we plot the disk mass as a function of a in Fig. [9] 
The EOS (and, equivalently, C) is the same for each plot. 
Again, we find a monotonic and steep increase of the disk 
mass as the increase of a for the fixed EOS and mass ra- 
tio. The enhancement of the disk mass by a progradc spin 
is more dramatic for the compact NS (for the soft EOS). 
For example, the difference in the disk mass between the 
cases of a = 0.75 and —0.5 is only by a factor of ~ 3 
when Q = 2 and 2H EOS is adopted. This low amplifi- 
cation is natural because tidal disruption of a large NS 
occurs at an orbit far enough from the ISCO for a sub- 
stantial amount of the disrupted material to escape the 
capture by the BH irrespective of a and because at such 
a large orbital separation the spin-orbit coupling effect 
is relatively weak. On the other hand, a few-orders-of- 



magnitude amplification of the disk mass is seen when 
Q = 2 and HB EOS is adopted. 

Finally, we comment on a possible unbound outflow. 
To estimate the rest mass of unbound material, we com- 
pute 



= / P* H {- 

J r>r A n 



-u t — l)d 3 x 7 



(34) 



where H(x) is a step function. Here, the material with 
Ut < — 1 should be considered to have an unbound orbit. 
We find that Af u b can be larger than 0.01M Q at ps 10 
ms after the merger for the stiff EOS like 2H and H, and 
a > 0. However, M u b does not approach a constant value 
and rather continues to decrease. Therefore, it is unclear 
whether M u b estimated at 10 ms after the merger can 
really become unbound or not, and we do not show the 
precise values of M u b- When the EOS is not stiff, M u b is 
negligible within the accuracy of our simulations. 



C. Structure of the remnant disk 

The structure of the remnant disk and its time evo- 
lution process depend on the mass ratio of the binary. 
We plot the rest-mass density profile at ps 5 and 10 ms 
after the onset of the merger for binaries with a = 0.75, 
HB EOS, and different values of Q in Fig. [TDJ The left 
column of Fig. [10] is plotted for ps 5 ms after the on- 
set of the merger and shows that the dense material of 
P ~ 10 9 s/ cm - 3 always extends to > 400 km. The spiral 
arm always spreads to a far region irrespective of EOSs, 
as far as the tidal disruption results in a massive disk. 
These plots also suggest that the accretion disk for a 
large value of Q — say, Q = 5 — keeps a nonaxisymmctric 
structure in the vicinity of the remnant BH at this time. 
This feature is qualitatively the same for binaries with 
other EOSs. When Q is small as ~ 2, the accretion disk 
becomes nearly axisymmetric in ps 5 ms after tidal dis- 
ruption because the dynamical time scale of the system 
(which is proportional to the BH mass) is shorter for a 
smaller value of Q. Also, because the ISCO radius of the 
BH is smaller, the maximum rest-mass density, p max , of 
the disk (which should be approximately proportional to 
the inverse square of the BH mass) reaches a higher value 
on average in time for a smaller value of Q. It should be 
noted that the difference in p max comes primarily from 
the difference in the radius and not from the difference 
in the disk mass, which do not vary by an order of mag- 
nitude for a = 0.75 and Q — 2-5. This difference in the 
nonaxisymmctric structure results in different features of 
gravitational waves (see Sec. IIV E[) . 

The middle and right columns of Fig. [TU] plot snap- 
shots at ps 10 ms after the onset of the merger. At this 
time, nonaxisymmetric structures are not as significant 
as those at ps 5 ms after the onset of the merger because 
the accretion disk settles toward an approximately sta- 
tionary state in the vicinity of the BH. The maximum 
values of the rest-mass density, /o max , in the accretion 
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TABLE V. Several key quantities for the merger remnants for Mns = 1.35Mq. All the quantities are estimated at ~ 10 ms 
after the approximate merger time t = imergcr. M r>rAH is the rest mass of the disk surrounding the BH; because the accretion 
is still ongoing due to the hydrodynamic angular momentum transport process, the listed values only give approximate masses 
of the long-lived accretion disks, which survive for a time longer than the dynamical time scale ~ a few ms. C e and C v are the 
circumferential radii of the AH along the equatorial and meridional planes, respectively. C e /4n and Mbh, c denote approximate 
masses of the remnant BH. Mj rr ,f is the irreducible mass of the remnant BH, and at is the nondimensional spin parameter of the 
remnant BH estimated from C p /C e . an and an are also the nondimensional spin parameters, estimated from the quantities on 
the AH and approximate conservation laws, respectively. We note that the values associated with the remnant BH for model 
B-Q2M135a75 (with an asterisk) are evaluated at ~ 5 ms after the onset of the merger because the BH area decreases by > 1% 
at « 10 ms after the onset of the merger and the error becomes large. 



Model 


M r>rAH [M G ] C e /47rM Mbh,o/Mq M irr , f /M C p /C e a f a f2 a n 


2H-Q2M135a75 
1.5H-Q2M135a75 
H-Q2M135a75 
HB-Q2M135a75 
B-Q2M135a75 


0.32 0.913 0.915 0.789 0.807 0.87 0.87 0.95 
0.29 0.918 0.920 0.785 0.794 0.89 0.89 0.95 
0.24 0.927 0.929 0.783 0.780 0.91 0.90 0.94 
0.21 0.933 0.934 0.783 0.772 0.91 0.91 0.94 
0.18 0.937* 0.938 0.790* 0.778* 0.91* 0.91* 0.93 


2H-Q2M135a5 
1.5H-Q2M135a5 
H-Q2M135a5 
HB-Q2M135a5 
B-Q2M135a5 


0.27 0.925 0.926 0.825 0.843 0.81 0.81 0.84 
0.23 0.935 0.936 0.831 0.840 0.82 0.81 0.84 
0.17 0.945 0.946 0.837 0.836 0.82 0.82 0.84 
0.14 0.951 0.952 0.840 0.832 0.83 0.83 0.84 
0.095 0.959 0.960 0.846 0.830 0.83 0.83 0.84 


2H-Q2M135a-5 
H-Q2M135a-5 

HB-Q2M135a-5 
B-Q2M135a-5 


0.13 0.961 0.962 0.931 0.954 0.48 0.48 0.50 
0.010 0.985 0.986 0.950 0.948 0.51 0.51 0.52 
0.0021 0.985 0.986 0.952 0.950 0.50 0.50 0.51 
2 x 10~ 4 0.983 0.984 0.952 0.952 0.49 0.49 0.50 


2H-Q3M135a75 
1.5H-Q3M135a75 
H-Q3M135a75 
HB-Q3M135a75 
B-Q3M135a75 


0.35 0.927 0.927 0.807 0.815 0.86 0.86 0.90 
0.30 0.931 0.934 0.811 0.815 0.86 0.86 0.90 
0.24 0.939 0.943 0.820 0.818 0.85 0.85 0.91 
0.22 0.941 0.943 0.812 0.805 0.87 0.87 0.90 
0.15 0.949 0.951 0.824 0.812 0.86 0.86 0.89 


2H-Q3M135a5 
1.5H-Q3M135a5 
H-Q3M135a5 
HB-Q3M135a5 
B-Q3M135a5 


0.28 0.939 0.940 0.858 0.874 0.74 0.74 0.77 
0.23 0.946 0.948 0.862 0.871 0.75 0.75 0.78 
0.16 0.955 0.957 0.867 0.866 0.76 0.76 0.78 
0.11 0.961 0.963 0.871 0.864 0.77 0.77 0.78 
0.050 0.969 0.971 0.877 0.862 0.77 0.77 0.79 


HB-Q3M135a-5 


< 10~ 4 0.986 0.987 0.973 0.980 0.32 0.32 0.33 


2H-Q4M135a75 
H-Q4M135a75 

HB-Q4M135a75 
B-Q4M135a75 


0.36 0.937 0.938 0.825 0.828 0.84 0.84 0.87 
0.23 0.948 0.951 0.831 0.823 0.84 0.84 0.88 
0.18 0.953 0.956 0.833 0.821 0.85 0.85 0.88 
0.11 0.960 0.963 0.837 0.817 0.85 0.85 0.88 


2H-Q4M135a5 
H-Q4M135a5 

HB-Q4M135a5 
B-Q4M135a5 


0.28 0.950 0.951 0.879 0.891 0.70 0.70 0.72 
0.085 0.970 0.973 0.890 0.880 0.73 0.73 0.74 
0.024 0.976 0.979 0.894 0.878 0.74 0.74 0.75 
0.0034 0.978 0.980 0.896 0.878 0.74 0.74 0.75 


2H-Q5M135a75 
H-Q5M135a75 

HB-Q5M135a75 
B-Q5M135a75 


0.36 0.946 0.947 0.838 0.835 0.82 0.82 0.85 
0.17 0.960 0.963 0.844 0.827 0.84 0.84 0.86 
0.095 0.966 0.970 0.848 0.824 0.84 0.84 0.86 
0.031 0.972 0.975 0.851 0.821 0.85 0.85 0.87 
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FIG. 8. The same as Fig. but for a = 0.75 and a = 0.5 with Q = 2-5. 
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FIG. 9. The disk mass, M r >r AH , at w 10 ms after the onset of the merger as a function of the nondimensional spin parameter 
of the BH a for models with Mns = 1.35Mq. The left and right panels show the results for models with the 2H and HB EOSs, 
respectively. The results revised from the previous work [3l[ are plotted for a = 0. 
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disk are still higher for a smaller value of Q. Indeed, the 
right column of Fig. QJJ shows that smaller values of Q 
result in producing a wider region with p > 10 12 g/cm 3 . 
By contrast, the disk for Q = 5 does not have such a high- 
density region. The smaller density may be unfavorable 
to be the short-hard GRB model. 

The size of a region where p > 10 10 g/cm 3 coincides 
approximately among four models with different values 
of Q and is always ~ 100 km. Furthermore, the mid- 
dle column suggests that the region of p > 10 s g/cm 3 
extends to larger distances when Q is larger. We plot 
the radial distribution of p along x and y axes for these 
models in Fig. [IT] Note that low-density regions near the 
origin are inside the BH. These plots show that p max is 
systematically higher for the binary with a smaller value 
of Q. These also show that the location of the isoden- 
sity surface of p = 10 10 g/cm 3 approximately coincides 
among different values of Q. Taking these facts into ac- 
count, we conclude that a typical profile of p(r) is steeper 
for smaller values of Q in the vicinity of the BH. A region 
of > 100 km away from the BH, where the profile p(r) 
shows relatively shallow decrease and < 10 10 g/cm 3 , cor- 
responds to the tail component, as is seen in the middle 
column of Fig. QJJJ 



D. Properties of the remnant BH 

Table [V] shows that masses and nondimensional spin 
parameters of the remnant BHs depend weakly on the 
adopted EOSs. The mass of the remnant BH tends 
to become large as the EOS softens for fixed values of 
(Q, Mns, a) for the case in which tidal disruption of the 
NS occurs. The reason for this is that the tidal disrup- 
tion occurs near the ISCO, and then the BH swallows a 
large amount of the NS mass when the EOS is soft. Ex- 
ceptionally, the mass of the BH becomes slightly larger 
for H and HB EOSs than for B EOS for binaries with 
(Q,a) = (2,-0.5). The reason for this is that the rem- 
nant disk masses are small as < O.OIMq for these cases 
and the amount of the energy radiated by gravitational 
waves primarily determines the final state (for more com- 
pact NSs, the radiated energy is larger because a closer 
inspiral orbit is achieved). The spin angular momen- 
tum of the remnant BH iSeH.f shows similar behavior to 
that of the BH mass. The situation becomes compli- 
cated for a spin parameter of the remnant BH defined 
by df = SsH.f A-^bh f> tne competition between the mass 
and angular momentum losses from the system makes 
the dependence of the nondimensional spin parameter of 
the remnant BH on the EOS very weak. For comparison, 
Ofi and af2 defined in Sec. 1111 Bl are also shown in Ta- 
ble |Vj As is found in our previous work (3l| , at and at2 
agree with each other within the error of Aa = 0.003. By 
contrast, an does not agree well with the other two esti- 
mates, as is found in another previous work of ours [2^ |. 
particularly when the massive remnant disk is formed 
and/or the mass of the BH is small: the maximum er- 



ror is Aa ~ 0.08. Taking into account the fact that the 
agreement between C e /47r and MeH.f is always better 
than 0.5%, a possible reason for this discrepancy is that 
Eq. (|3ip systematically underestimates the angular mo- 
mentum of the disk. Hereafter, we only refer to at as the 
nondimensional spin parameter of the remnant BH. 

The nondimensional spin parameter of the remnant BH 
depends strongly on the initial spin parameter, a, and 
the mass ratio, Q. Approximate values of the nondimen- 
sional spin parameter of the remnant BH, at , are shown 
in Fig. [12] as a function of the initial BH spin parameter, 
a. We also plot lines obtained by a linear fitting of data 
for Q = 2 and 3 of the following form, 

at = 0.32a + 0.66 (Q = 2), (35) 
Of = 0.43a + 0.54 (Q = 3). (36) 

The relation for Q = 3 agrees approximately with the 
results reported in Refs. (28|,[85| within an error of < 5%, 
and the agreement becomes better for a larger value of a. 
Figure [T^] and these relations show that af is an approx- 
imately linear function of a. In a zeroth approximation, 
the slope and intercept of the linear relation denote the 
contribution from the initial BH spin angular momen- 
tum, Sbh, and the orbital angular momentum of the bi- 
nary, Jo , respectively. The larger slope for a larger value 
of Q is explained by a larger contribution from the spin 
of the initial BH of mass A/bh = Q-Mns to the spin of the 
remnant BH of mass MBH,f ~ (1 + Q)Mns- These pre- 
dict the value of the slope to be Q 2 / (1 + Q) 2 . However, 
the slope obtained by numerical simulations is smaller by 
~ 25-30 % than this predicted slope, because the amount 
of angular momenta redistributed to the remnant disk 
and extracted by gravitational waves become larger for a 
larger value of a. The fitting function also suggests that 
the merger of an extremely spinning BH of a = 1 and a 
NS with an irrotational velocity field results in a remnant 
BH with at sa 0.98 for BH-NS binaries with Q = 2 and 
3 and hence never forms an overspinning BH, i.e., a BH 
with af > 1. Furthermore, the results for Q = 4 shown 
in Table fVl also suggest af « 0.97 for the merger of an 
extremely spinning BH and an irrotational NS. These re- 
sults give a circumstantial support for cosmic censorship 
conjecture 86]. Whether af < 0.98(< 1) is an universal 
consequence of a general BH-NS binary merger or not 
should be confirmed by simulations of higher mass-ratio 
binary mergers, in particular, with (nearly) extremal BH 
spin. 

From these typical values of af and Mbh,£ j we can esti- 
mate typical frequencies of quasinormal modes (hereafter 
QNM) /qnm of the remnant BH by a fitting formula (8?| 

/QNMA/BH.f « ^-[1.5251 - 1.1568(1 - a) - 1292 ]. (37) 

We plot these values in Fig. [T3] They are in good agree- 
ment with those of the ringdown waveforms, if the QNM 
is excited after the merger. 
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FIG. 10. The same as Fig.[3]with contour curves for p = 10 10 and 10 12 g/ cm plotted. In all the plots, HB EOS and a = 0.75 are 
adopted. The first, second, third, and fourth rows are for Q = 2, 3, 4, and 5, respectively. The left column plots the snapshots 
at 5 ms after the onset of the merger. The middle column plots the snapshots at 10 ms after the onset of the merger, and the 
right column plots close-ups of the middle column. 
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FIG. 11. Radial distribution of the rest-mass density at 10 ms after the onset of the merger for different values of Q. The left 
and right panels show the distribution along the x and y axes, respectively. In both plots, (Mns,o) = (1.35M0 , 0.75) and HB 
EOS are adopted. 
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FIG. 12. The typical nondimensional spin parameters of the 
remnant BH at as a function of the initial BH spin parameter 
a. The solid lines are obtained by a linear fitting of the data 
for Q = 2 and Q = 3. 
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FIG. 13. The typical QNM frequency of the remnant BH 
normalized by its mass /cjNMA/BH.f. 



E. Gravitational waveforms 

In this section, we show (Z,m) = (2,2), plus-mode 
gravitational waveforms h+ for selected models obtained 
in this study, as well as the waveform for models obtained 
in our previous simulations 31 1. We plot all the wave- 



forms for an observer along the z axis as a function of the 
approximate retarded time 



t Iet =t-D-2MoH D / M o) 



(38) 



The amplitude of the waveforms is normalized as 
Dh + /m or we show physical amplitude observed at a 
hypothetical distance D = 100 Mpc along the z axis. 
Gravitational waveforms calculated in the Taylor-T4 for- 
mula are plotted together in the figures to validate the 
waveforms obtained in our numerical simulations during 
the inspiral phase. Numerical waveforms during 2-3 ini- 
tial cycles deviate from ones obtained from the Taylor- T4 
formula in all the cases due to the lack of an approaching 
velocity in the initial data. This deficit is ascribed to in- 
sufficient modeling of the quasiequilibrium state and im- 
provement in the future is important to obtain more ac- 
curate gravitational- wave templates [H, |89[ . Our wave- 
forms also deviate from the Taylor-T4 waveforms in the 
late inspiral phase due to a physical reason, which we de- 
scribe below. Comparisons between waveforms obtained 
from simulations with different grid resolutions are shown 
in the Appendix „ 

Figure Q3] shows the gravitational waveforms for bina- 
ries with HB EOS, (Q, M NS ) = (3,1.35M Q ) but with 
different BH spin parameters, a = 0.75, 0.5, 0, and 
—0.5. This figure shows that the time to the merger, 
to which we refer approximately as the time at which 
the maximum gravitational-wave amplitude is achieved, 
for n m = 0.030 becomes longer by ps 10 ms as the in- 
crease of the BH spin within the range concerned here. 
This difference in the merger time owes primarily to the 
spin-orbit interaction described in Sec. IIV A[ and this 



21 




o 

Q. 



O 
O 




O 
Q_ 



O 
O 




o 

Q. 



O 
O 




5e-22 
4e-22 
3e-22 
2e-22 
1e-22 


- -1e-22 
-2e-22 
-3e-22 
-4e-22 
-5e-22 



30 



o 

Q_ 



O 

o 



FIG. 14. (l,m) = (2,2), plus-mode gravitational waves for models HB-Q3M135a75, HB-Q3M135a5, HB-Q3M135, and HB- 
Q3M135a-5. All the waveforms are shown for an observer located along the z axis (the axis of rotation) and plotted as a 
function of a retarded time. The left axis denotes the amplitude normalized by the distance from the binary D and the total 
mass mo. The right axis denotes the physical amplitude of gravitational waves observed at a hypothetical distance 100 Mpc. 
The dotted curves denote the waveform calculated by the Taylor-T4 formula. 



behavior is qualitatively the same for binaries with any 
EOS. The numerical and Taylor- T4 waveforms agree well 
with each other during an inspiral phase for all the cases. 

For the progradc BH spin cases, the Taylor-T4 for- 
mula does not track the evolution for ~ 0.5 inspiral or- 
bit just before the merger. The Taylor- T4 amplitude 
departs from that of numerical relativity and even di- 
verges. Accordingly, the number of gravitational- wave 
cycles differs by as much as unity between the numer- 
ical and Taylor- T4 waveforms. The difference in the 
number of cycles is larger for higher mass-ratio binaries 
with progradc BH spins. We show the waveforms for 
binaries with (Q,A/ NS ,a) = (4, 1.35M©, 0.75) and with 
(Q,M NS ,a) = (5,1.35Af Q ,0.75) for 2H, H, HB, and B 
EOSs in Figs. [15] and [T51 respectively The deviation is 
clear for H, HB, and B EOSs in both figures. This dif- 
ference indicates that the phase evolution predicted by 
the Taylor-T4 formula is not sufficient to model the last 
inspiral phase of a coalescing binary with the high mass 
ratio of Q > 3 and the prograde BH spin of a > 0.5. 

For the retrograde BH spin case (the bottom right 
panel of Fig. [T3j) . the phase evolution deviates between 
the numerical and Taylor- T4 waveforms in the last orbit 



before the merger. This deviation may be partly ascribed 
to the small number of orbits in our simulation but ap- 
pears to be primarily ascribed to a larger angular veloc- 
ity, or equivalently a larger PN parameter, f2mo, at the 
last orbit for a retrograde BH spin. Thus, the Taylor- T4 
formula seems to be again insufficient for modeling the 
retrograde BH spin cases. 

Figure [14] also shows that the gravitational waveform 
in the merger stage depends strongly on the BH spin. 
For a binary with (Q,a) = (3,0.75), gravitational waves 
show a sudden decrease in the amplitude at t re t ~ 27 ms, 
which is a clear signature of tidal disruption. Gravita- 
tional waves associated with the ringdown of a remnant 
BH are absent due to the phase cancellation by nearly 
axisymmetric accretion of the disrupted material. This 
feature is consistent with the formation of a massive rem- 
nant disk, which is described in Sec. IIV B[ for the pro- 
grade BH spin. For binaries with (Q, a) = (3, < 0), on the 
other hand, gravitational waves end up with ringdown 
waveforms associated with the remnant BHs because the 
tidal effect is very weak throughout the merger. In these 
circumstances, gravitational waves do not show strong 
signatures of tidal deformation and disruption of the NS. 
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FIG. 15. The same as Fig. [HI but for models 2H-Q4M135a75, H-Q4M135a75, HB-Q4M135a75, and B-Q4M135a75. 



Gravitational waves for a binary with (Q, a) = (3, 0.5) 
show a qualitatively new feature (the top right panel of 
Fig. 1141) . In this case, a ringdown waveform of the rem- 
nant BH is seen in the final stage, although the NS is 
tidally disrupted and the disk mass is larger than O.IMq. 
Namely, both tidal disruption of the NS and excitation 
of a QNM of the remnant BH occur in a compatible 
manner. The same feature is also found for a binary 
with a high mass ratio and a prograde BH spin, i.e., 
(Q,a) = (> 4,0.75), shown in Figs. [TBI and [TBI except 
for 2H EOS, with which the NS is disrupted at a fairly 
distant orbit. These waveforms are often seen for BH-NS 
binaries with a heavy BH (or a high mass ratio) with the 
prograde BH spin, which results in the NS tidal disrup- 
tion, and is never seen for BH-NS binaries with Q = 2 or 
high mass-ratio binaries with nonspinning BHs. 

The ratio of the areal radius of the remnant BH to 
the NS radius, Rns, is intimately related to the differ- 
ent excitation degree of the QNM between low and high 
mass-ratio binaries in the presence of NS tidal disrup- 
tion (9?j, . Schematic pictures of merger processes are 
depicted in Fig.[T7J If tidal disruption does not occur, the 
NS is simply swallowed by the BH and excites a QNM, as 
shown in the middle panel of Fig. 1171 If tidal disruption 
occurs in a binary with a low mass ratio, the disrupted 
material spreads around the BH to soon form a nearly 



axisymmctric disk. Approximately speaking, this occurs 
if the BH areal radius is smaller than the NS radius, 
as is shown in the left panel of Fig. [T7J Thus, the NS 
tidal disruption has a strong effect to suppress the exci- 
tation of a QNM through the phase cancellation in the 
low mass-ratio binary. However, the situation is differ- 
ent in a high mass-ratio binary. Whereas the disrupted 
material forms an axisymmetric accretion disk around 
the BH in a sufficiently long time duration, the accretion 
just after the merger does not proceed in an axisymmet- 
ric way in high mass-ratio binaries, such as Q = 4, except 
for the extremely stiff EOS. This is because the BH ra- 
dius for Q = 4 approximately doubles that for Q = 2, and 
hence, the disrupted material takes longer time to spread 
around the BH. In other words, the NS material accretes 
onto the BH coherently even after the tidal disruption, 
as is shown in the right panel of Fig. [TTl because the 
BH radius is so large that the disrupted material can- 
not fully cover the BH surface before the BH swallows 
a large portion of the material. In the exceptional 2H 
EOS case, tidal disruption occurs sufficiently far outside 
the BH due to the large radius of the NS, and hence, 
the disrupted material is able to spread around the BH 
to form a nearly axisymmetric accretion disk before the 
prompt infall. Therefore, the QNM of a remnant BH is 
not excited for 2H EOS. 
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FIG. 16. The same as Fig. [HI but for models 2H-Q5M135a75, H-Q5M135a75, HB-Q5M135a75, and B-Q5M135a75. 




FIG. 17. Schematic pictures for three types of the merger process. The solid filled circle denotes the BH, the distorted ellipsoid 
denotes the NS, the solid circle is the location of the ISCO, and the dashed circle is the location of the radius at which the 
tidal disruption occurs. Left: the NS is tidally disrupted, and the spatial extent of the disrupted material is larger than or as 
large as that of the BH. Middle: the NS is not tidally disrupted. Right: the NS is tidally disrupted, and the spatial extent of 
the disrupted material is smaller than that of the BH. 



F. Gravitational- wave spectrum spectrum as a sum of each Fourier component of two 



Key features of gravitational waves are reflected in the 
Fourier spectrum. In this paper, we define the Fourier 
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independent polarizations of the (I, |m|) = (2, 2) mode as 




M/)I 2 + IM/)I 2 



27lift h A (t)dt, 



(39) 
(40) 



where A denotes two polarization modes, + or x . We 
show a nondimcnsional spectrum, fh(f), observed at a 
hypothetical distance of 100 Mpc as a function of the 
gravitational- wave frequency, / (Hz), or a normalized 
amplitude, fh(f)D/mo, as a function of a nondimcn- 
sional frequency, /mo. The amplitude of gravitational 
waves, h,Ai is given as the amplitude observed along the 
z axis, which is the most optimistic direction for the 
gravitational-wave detection. We note that the actual 
amplitude of gravitational waves depends on an angle lo- 
cating the source in the sky and on an angle specifying 
the orientation of the orbital plane of the binary. The 
angular average of the effective amplitude is w OAfh(f). 
We always exclude spurious radiation components for 
i rc t < ms ) using a step function of the retarded time as 
a window function. 

To show the dependence of the gravitational-wave 
spectra on the BH spin parameter, we plot the spectra for 
models HB-Q2M135a75, HB-Q2M135a5, HB-Q2M135, 
and HB-Q2M135a-5 in the left panel of Fig. [TJ and for 
models HB-Q3M135a75, HB-Q3M135a5, HB-Q3M135, 
and HB-Q3M135a-5 in the right panel of Fig. [U In 
the early inspiral phase of / < 1 kHz, where the point- 
particle approximation works well, the amplitude of the 
gravitational-wave spectrum for a given frequency in- 
creases monotonically as a increases. This is a feature 
expected from the PN calculation and is explained by the 
spin-orbit interaction as follows: The power spectrum of 
gravitational radiation is written as 
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(41) 



On the other hand, retaining only 1.5PN, the lowest- 
order spin-orbit interaction terms, Eqs. (4.10) and (4.14) 
of Ref . [79| , derive the expression for this quantity as 
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(42) 

Thus, the effective amplitude, fh(f), for a given fre- 
quency / increases monotonically as the BH spin param- 
eter, a, increases in the inspiral phase. 

Figure [U plots the spectra for models 2H-Q2M135a75, 
H-Q2M135a75, HB-Q2M135a75, and B-Q2M135a75, for 
which only the EOS is different, and indicates that the 
amplitude in the early inspiral phase does not depend 
strongly on C. This is because the finite-size effect of the 
NS does not play an important role in the early inspi- 
ral phase (but, see Refs. [93l . |94|0 . as already found for 
nonspinning BH-NS binaries [31] . 



In the late inspiral phase of 1 kHz < / < / cut , where 
/ cut is a characteristic frequency at which the spectrum 
starts damping exponentially (see below), the amplitude 
is significantly larger than the Taylor- T4 formula for the 
cases in which the NS is not disrupted. This is because 
the binaries in the inspiral and plunge after the NS en- 
ters the BH's ISCO emit gravitational waves in reality, 
whereas the Taylor- T4 formula docs not take into ac- 
count the motion inside the ISCO. In contrast to the 
spectrum calculated by the Taylor-T4 formula, which de- 
creases steeply after the last inspiral phase, the amplitude 
obtained from the simulation depends only weakly on the 
gravitational- wave frequency in that phase, as far as the 
tidal disruption does not occur. 

The most fruitful information of the NS comes from 
the gravitational-wave spectrum in the merger phase 
through the "cutoff frequency," / cut , which depends on 
the BH spin as well as the NS compactness [2!| HH . If 
the NS tidal disruption occurs, the spectrum damps at 
/ = /tidal ~ 2-4 kHz, which denotes the frequency at 
the tidal disruption and depends sensitively on physi- 
cal parameters of the binary. In that case, gravitational 
waves for a higher frequency, / > / cut ~ /tidal, are not 
emitted by the binary in the inspiral motion but only 
weakly by disrupted material. Because the disrupted 
material gradually spreads around the BH to form a 
nearly axisymmetric disk, the emission of gravitational 
waves is suppressed at the high frequency. Thus, the 
spectrum shows a relatively moderate damping around 
/ ~ /cut, which is closely related to the NS compactness 
through the tidal disruption. The spectra for binaries 
with (Q, a) = (2, > 0) and (3, > 0.5) in the left panel 
of Fig. [15] correspond to these cases. We see that the 
cutoff frequency, / cu t, for these models decreases as the 
BH spin parameter increases. This is ascribed to the de- 
crease of the orbital frequency at the tidal disruption for 
a binary with the prograde BH spin. The enhancement 
of the effective centrifugal force by the spin-orbit inter- 
action reduces the orbital frequency at the tidal disrup- 
tion, /tidal, although the orbital separation at the tidal 
disruption itself does not vary much even in the presence 
of the BH spin. If the tidal disruption does not occur 
during the merger, however, inspiral-like motion contin- 
ues at higher frequencies near and even inside the ISCO 
until the BH swallows the NS. In this case, the spectrum 
amplitude depends only weakly on / in the frequency 
range / < / cu t and damps for / > / cut , which is closely 
related to the QNM frequency of the remnant BH, /qnm- 
The spectra for (Q, a) = (2, -0.5) and (3, < 0) in Fig. [IS] 
show this feature. Note that the amplitude for model 
HB-Q3M135a-5 is smaller than for model HB-Q3M135 
for the frequency range shown in Fig. [TS] because tidal 
disruption does not play an important role, and Eq. (|42|) 
applies throughout the merger in both cases. 

It is noteworthy that a prograde BH spin is favorable 
for the gravitational- wave detection in the inspiral phase 
and the estimation of / cu t in the merger phase because 
the prograde spin enhances the amplitude for a given 
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FIG. 18. Gravitational-wave spectra for BH-NS binaries with HB EOS, Af N s = 1.35M and a = 0.75, 0.5, 0, and -0.5. The 
left and right panels show the spectra for Q = 2 and 3, respectively. The upper axis denotes the normalized frequency, /mo, 
and the right axis denotes the normalized amplitude, fh(f)D/mo. The bottom axis denotes the frequency, /, in Hz, and the 
left axis denotes the nondimensional amplitude of gravitational waves, fh(f), observed at a hypothetical distance 100 Mpc 
from the binary along the z axis. The dashed curves are planned noise curves of the LCGT ("LCGT"), the Advanced LIGO 
optimized for IAMq NS-NS detection ("Standard"), the Advanced LIGO optimized for the burst detection ("Broadband"), 
and the Einstein Telescope ("ET") [H . 
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FIG. 19. The same as Fig. [TS] but for (Q,M NS ,a) = 
(2, 1.35M Q ,0.75) with 2H, H, HB, and B EOSs. The spec- 
trum derived by the Taylor- T4 formula is also included. 



ries (except for the model with 2H EOS), both the NS 
tidal disruption and the excitation of a QNM of the rem- 
nant BH occur as is described in Sec. IIVEI Hence, the 
gravitational spectrum has two characteristic frequencies, 
i.e., /tidal and /qnm, simultaneously. The spectra plot- 
ted in Fig. [20] indeed show such features. After the NS 
is tidally disrupted, the amplitude of the gravitational- 
wave spectrum shows a slow damp for / > /tidal ~ 2 kHz. 
Then, the spectrum damps steeply above the frequency 
of the QNM, / > /qnm ~ 3 kHz. A schematic figure of 
different spectra is depicted in Fig. [2U and the spectrum 
described in this paragraph corresponds to spectrum (iii) 
in this figure. This suggests that the cutoff frequency, 
/cut, of a high mass-ratio binary is not determined by 
a unique physical process like NS tidal disruption or a 
ringdown of a remnant BH, as far as both of them occur. 

To estimate the cutoff frequency quantitatively, we fit 
the gravitational-wave spectra by a function with seven 
free parameters of the form 



frequency in the inspiral phase and decreases the cutoff 
frequency in the merger phase. Note that the most sensi- 
tive frequency range for ground-based detectors is / ~ 10 
Hz-1 kHz, which is usually lower than / cu t- Thus, the 
features found here are encouraging for the gravitational- 
wave astronomy to become an important tool for inves- 
tigating the NS radius and EOS. 

The gravitational-wave spectra of binaries with high 
mass ratios show qualitatively different behavior for 
a high frequency. Figure [501 plots the gravitational- 
wave spectra obtained for models with (Q, M^s,a) = 
(4, 1.35M , 0.75) and (5, 1.35M Q , 0.75). For these bina- 
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where /i3pn(/) is the Fourier spectrum calculated by the 
Taylor-T4 formula. The first term in Eq. (j43|) models 
the inspiral spectrum, and the second term models the 
merger and ringdown spectra. Wc determine seven free 
parameters / ins , / ins2 , ./dam, 0" ins , cr ins2 , a dam> and A by 
the condition that the following weighted norm is mini- 
mized: 
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FIG. 21. A schematic figure of three types of gravitational- 
wave spectra. Spectrum (i) is for the case in which tidal dis- 
ruption occurs far outside the ISCO, and spectrum (ii) is for 
the case in which tidal disruption does not occur. Spectrum 
(iii) is for the case in which tidal disruption occurs and the 
QNM is also excited. The filled and open circles denote /tidal 
and /qnm, respectively. 



Here, i denotes the data point for the spectrum. In the 
previous works [2!| HH , we identify /dam in Eq. with 
/ cut , which is most strongly correlated with the NS com- 
pactness for nonspinning BH-NS binaries [9j|. In the 
present work, however, we obtain no strong correlation 
between /dam (and the other parameters) and any param- 
eter of physical importance, such as a or C. The reason 
may be ascribed to the inadequacy of the functional form 
of Eq. ([4"5]l . where the set of free parameters is degener- 
ate to some extent. In particular, such a degeneracy is 
severe for a high mass-ratio binary due to two reasons. 
First, modeling an inspiral spectrum by the Taylor-T4 
formula is inadequate for the late inspiral phase of a high 
mass-ratio binary due to the lack of information from the 
Taylor- T4 formula, as is described in Sec. IIV El Second, 



FIG. 22. The fitting for model HB-Q2M135a75. The long- 
dashed and middle-dashed curves show the first and second 
terms of Eq. (|43|l , respectively. The short-dashed curve is the 
sum of these two terms. 



there is no unique, physically motivated identification of 
/cut when both the NS tidal disruption and the QNM ex- 
citation occur. (Fortunately, these degeneracies did not 
cause problems in the case of the nonspinning BH-NS 
binary with a low mass ratio [3l|.) To overcome these 
problems with the fitting procedure, we redefine / cu t as 
the higher one of two frequencies at which the second 
term in Eq. (|43|) takes a half value of its maximum. An 
example of this fitting procedure is shown in Fig. 1221 In 
this figure, ff max corresponds to the maximum value of 
the second term in Eq. (|43p . We find that this definition 
of /cut works well to read off the NS compactness from 
the gravitational-wave spectrum. 

Figure [2^1 shows /cut^o for spectra obtained for all bi- 
naries with Q = 2 as a function of the NS compactness, 
C, in logarithmic scales. We also plot the typical QNM 
frequency of the remnant BH, /qnm, which depends pri- 
marily on a of the initial BH for a fixed value of Q. For 
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FIG. 23. The cutoff frequency times the total mass / cu tmo 
as a function of the NS compactness C for Q — 2 binaries in 
logarithmic scales. The solid and dashed lines are obtained 
by linear fittings of data for a — 0.75 and a = 0, respectively. 
Horizontal lines denote the typical QNM frequencies of the 
remnant BHs. We note that /qnm > /cut for a — 0.75 as long 
as C < 0.2. 



each value of a, we find that /cut^o increases monotoni- 
cally as C increases, and an approximate power law holds 
as 

ln(/ cut m ) »p(fl) hi C + q(a) (Q = 2), (45) 

where p(a)(> 0) and q(a) depend only on a, for any 
value of a when Q = 2. This monotonic relation be- 
tween /cut^o and C suggests us a possibility to extract 
the compactness, C, of a NS from the gravitational- wave 
observation. It is noteworthy that this relation includes 
only C but neither A/ns nor -Rns independently. It should 
also be noted that the simple relation found here is a 
consequence of our choice for a common value of the adi- 
abatic index of the core EOS, T2 = 3 [3l| . The increase 
of /cut"io with the increase of C indicates that a more 
compact NS is less subject to the BH tidal effect and dis- 
rupted at a closer orbit to the BH than a less compact 
NS is. The difference in /cut^o due to the difference in 
a becomes clearer for larger values of C, and conversely, 
/ cut mo depends only weakly on a if the compactness is 
as small as ~ 0.12. The weak dependence on a for the 
small values of C is due to the fact that the effect of the 
BH spin at a distant orbit, at which the NS with a large 
radius is disrupted, is weak. 

Figure [23] also shows that p(a) is a decreasing function 
of a. More specifically we obtain the relations 

ln(/cutm ) = (2.92 ± 0.06) InC + (2.32 ± 0.12) (46) 

for a = [96[ and 

ln(/ cut mo) = (2.39 ± 0.06) InC + (1.11 ± 0.11) (47) 

for a — 0.75 by a linear fitting. The decreasing nature of 
p(a) is explained by the fact that the spin-orbit repulsive 



force for the progradc BH spin, which reduces the orbital 
frequency at the NS tidal disruption, works efficiently for 
a close orbit and hence for the NS with a large value of 
C. It is important that p(a) is always larger than 1.5, 
which is expected from the analysis of the condition for 
the mass shedding, 



C 3 / 2 (l + Q) 3 / 2 



(48) 



The large value of p(a) is favorable for determining the 
NS compactness from the gravitational- wave observation 
because the dependence of / cu t on C becomes stronger. 
Note that / cu t is always lower than the QNM frequency of 
the remnant BH for a realistic range of the compactness 
C < 0.2 for (Q, a) = (2, > 0). If a is negative, on the other 
hand, / cu t for the binary with a compact NS of C > 0.18 
may be determined by the QNM frequency, /qnm, and 
it will be difficult to determine the NS compactness from 
the cutoff frequency. 

Figure [Ml shows the f C ntmo~C relation of gravitational- 
wave spectra obtained for all binaries with a = 0.75 and 
a = 0.5. This figure, combined with Fig. [23l clearly 
indicates that the approximate power law of the form 



ln(/ cut m ) = p(Q, a) InC + q(Q, a) 



(49) 



holds for binaries of C < 0.2 with Q = 5 as far as a ~ 0.75 
and with Q < 4 as far as a ~ 0.5. The striking feature 
is that the cutoff frequency is lower than the QNM fre- 
quency of the remnant BH, /qnm, for (Q, a) = (< 4, 0.75) 
and for (Q, a) = (< 3, 0.5) even if a QNM is excited. Ac- 
cordingly, / cut shows a strong correlation with C. For 
(Q,a) = (5,0.75) and (4,0.5), / cut is lower than /qnm as 
far as C < 0.18 and 0.17, respectively, and, therefore, the 
strong correlation between / cut and C is found within this 
range. Although / cut for the binary with a high mass ra- 
tio should not be considered as /tidal due to the QNM ex- 
citation, monotonic relations between / cu t?7io and C gives 
us an opportunity to explore the NS radius and EOS. It 
should be noted that gravitational waves from a higher 
mass-ratio binary are more subject to the gravitational- 
wave detection due to the larger amplitude in the inspiral 
phase and the lower cutoff frequency. We again note that 
a massive BH of Mbh ^ 5M© is an astrophysically real- 
istic consequence of the stellar evolution [H, H3] • Taking 
these facts into account, we conclude that gravitational 
waves from the BH-NS binary are a promising tool to 
investigate the NS radius and EOS in the next decade. 



G. Energy and angular momentum radiated by 
gravitational waves 

Table fVTl lists the total energy AE/Mq and angular 
momentum A J / Jq radiated by gravitational waves. We 
estimate systematic errors in the estimation of AE and 
A J to be ~ 10%, which arc ascribed mainly to the finite 
grid resolution and partly to the finite extraction radius 



28 



a=0.75 



a=0.5 




0.12 0.14 0.16 0.18 0.2 
Compactness 




0.12 0.14 0.16 0.18 0.2 
Compactness 



FIG. 24. The same as Fig. [23]but for a = 0.75 (left) and 0.5 (right). In both panels, /qnm > /cut for Q = 2. 



for 4"4. Because AE and A J depend on the choice of 
Uomoj we do not compare directly the results obtained 
for models with different values of i~tomo and accordingly 
models with different values of Q (see Tabic ITTT|) . 

Contributions from all the I = 2-4 modes of gravita- 
tional waves are taken into account. The (I, \m\) = (2, 2) 
mode always contributes by > 85% to AE and AJ. AE 
and A J taken away by higher-mode gravitational waves 
are substantial for high mass-ratio binaries. For example, 
the (3, 3) mode contributes by ~ 2, 5, 7.5, and 10% for 
binaries with Q = 2, 3, 4, and 5, respectively. The (4, 4) 
mode gravitational waves contribute by 1 ~ 2% for bina- 
ries with Q = 3-5. These values depend only weakly on 
a and the EOS, and contributions of modes with I =/= m 
arc negligible compared to those of I = m modes. 

Table [VI] shows that AE/M and AJ/J increase 
monotonically as the NS compactness, C, increases for 
binaries with fixed values of (Q, a). This is the same re- 
sult as that obtained for nonspinning BH-NS binaries [3l[ 
and is explained by a longer inspiral phase for a softer 
EOS due to the later onset of mass shedding and the 
later tidal disruption. The ratio between these two val- 
ues, (AJ/Jq)/(AE/Mq), decreases as C increases. This 
agrees again with the result for the nonspinning BH cases 
and is explained by a relation AJ / AE w m/fi for a given 
angular harmonic of m and by the fact that more radia- 
tion is emitted from the orbit of a larger value of f2 for 
a softer EOS. Note that these arguments are based on 
little dependence of gravitational-wave luminosity in the 
inspiral phase on C for a fixed value of a; tidal correction 
to the luminosity in the inspiral phase is not important. 

Table IVT1 shows that AE/Mq does not depend strongly 
on a, while AJ/ Jo increases as a increases in many cases 
for a fixed value of C. Remember that dE/df in the in- 
spiral phase increases for a large value of a, as is given by 
Eq. (|42j) . However, the orbital frequency, f2, at the tidal 
disruption decreases for a large value of a due to the spin- 
orbit interaction. Because of these two competing effects, 
the binding energy at the tidal disruption depends only 



weakly on a, and hence, AE/Mq does not change very 
much among different values of a. The increase of AJ/ Jo 
for a large value of a is due to the large value of dE/df 
in the inspiral phase, during which Q is relatively low, 
and to the approximate relation AJ « mAE/Q, which 
enhances the contribution of low-frequency gravitational 
waves. 

Finally, we comment on the linear momentum AP ra- 
diated by gravitational waves and an associated kick ve- 
locity fkick = AP/Mo of the remnant BH. Because of the 
mass and spin asymmetries, the remnant BH achieves 
the kick velocity of ~ 100-250 km/s when the effect of 
tidal disruption is weak, e.g., for models HB-Q3M135a- 
5 and B-Q3M135. Although our results for AP do not 
converge as well as those for AE and A J due to the slow 
convergence of (l,m) ^ (2,2) mode gravitational waves, 
the values of ^kick are in reasonable agreements with the 
fitting formula derived using the results of simulations 
for the binary BH merger [93, H3|. By contrast, fkick is 
suppressed to < 100 km/s when tidal disruption occurs 
far outside the ISCO. The reason for this is that the tidal 
disruption suppresses significantly the gravitational radi- 
ation from the last inspiral and merger phases, during 
which the linear momentum is emitted most efficiently 
This trend is consistent with the result found in our pre- 
vious work [29l ]. 



V. SUMMARY 

We performed numerical simulations for the merger 
of BH-NS binaries with various BH spins aligned or an- 
tialigned with the orbital angular momentum, using an 
AMR code SACRA with systematically chosen five piece- 
wise polytropic EOSs. We investigated the dependence of 
the merger process, properties and structures of the rem- 
nant disk, properties of the remnant BH, gravitational 
waveforms, and their spectra on the spin of the BH and 
the EOS of the NS. In particular, we focused on the case 
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TABLE VI. Total radiated energy AE and angular momen- 
tum AJ carried away by gravitational waves. AE and AJ 
are normalized with respect to the initial ADM mass Mo and 
angular momentum Jo, respectively. We also show the ratio 
between AJ and AE. 



Model 


AE/M 


(%) AJ/Jo 


(%) (AJ/J )/{AE/M ) 


2H-Q2M135a75 


0.58 


16 


27 


1.5H-Q2M135a75 


0.79 


19 


24 


H-Q2M135a75 


1.1 


24 


21 


HB-Q2M135a75 


1.4 


26 


19 


B-Q2M135a75 


1.7 


29 


17 


2H-Q2M135a5 


0.60 


17 


26 


1.5H-Q2M135a5 


0.79 


19 


24 


H-Q2M135a5 


1.2 


24 


20 


HB-Q2M135a5 


1.4 


26 


19 


B-Q2M135a5 


1.7 


28 


17 


2H-Q2M135a-5 


0.57 


15 


26 


H-Q2M135a-5 


1.1 


19 


16 


HB-Q2M135a-5 


1.4 


19 


14 


B-Q2M135a-5 


1.6 


21 


13 


2H-Q2M12a75 


0.40 


12 


30 


H-Q2M12a75 


0.79 


19 


24 


HB-Q2M12a75 


0.95 


21 


22 


B-Q2M12a75 


1.2 


24 


21 


2H-Q2M145a75 


0.73 


19 


25 


H-Q2M145a75 


1.5 


27 


19 


HB-Q2M145a75 


1.7 


30 


17 


B-Q2M145a75 


2.1 


32 


15 


2H-Q3M135a75 


0.72 


20 


28 


1.5H-Q3M135a75 


0.97 


23 


24 


H-Q3M135a75 


1.3 


27 


20 


HB-Q3M135a75 


1.6 


30 


19 


B-Q3M135a75 


2.0 


34 


17 


2H-Q3M135a5 


0.70 


19 


27 


1.5H-Q3M135a5 


0.94 


22 


23 


H-Q3M135a5 


1.4 


26 


19 


HB-Q3M135a5 


1.7 


29 


17 


B-Q3M135a5 


2.0 


31 


15 


HB-Q3M135a-5 


1.3 


19 


14 


2H-Q3M145a75 


0.88 


22 


25 


H-Q3M145a75 


1.7 


31 


18 


HB-Q3M145a75 


2.1 


34 


16 


B-Q3M145a75 


2.5 


37 


15 


2H-Q4M135a75 


0.81 


23 


28 


H-Q4M135a75 


1.5 


31 


21 


HB-Q4M135a75 


1.8 


33 


19 


B-Q4M135a75 


2.1 


36 


17 


2H-Q4M135a5 


0.72 


19 


27 


H-Q4M135a5 


1.5 


27 


19 


HB-Q4M135a5 


1.7 


29 


17 


B-Q4M135a5 


1.9 


31 


16 


2H-Q5M135a75 


0.83 


24 


29 


H-Q5M135a75 


1.6 


33 


20 


HB-Q5M135a75 


1.!) 


35 


19 


B-Q5M135a75 


2.1 


36 


17 



in which the BH has a progradc spin, and the tidal dis- 
ruption of the NS by a companion BH plays an important 
role. We adopted a number of parameters for the mass 
ratio, NS mass, and BH spin. By preparing the initial 
condition with a distant orbit and a small eccentricity, 
we always tracked > 5 quasicircular orbits in the inspiral 
phase and studied the merger phase in a realistic setting. 
The treatment of hydrodynamic equations and the es- 
timation method of the disk mass are improved in this 
work. In the following, we summarize the conclusions in 
this paper: 

1. It is shown that a progradc BH spin enhances the 
effect of NS tidal disruption by the spin-orbit in- 
teraction. The mass of the remnant disk increases 
as the BH spin increases because the ISCO radius 
of the BH becomes small. A remarkable point is 
that the BH-NS binary with a high mass ratio of 
even Q = 5 can form a sufficiently massive disk 
of > 0.1 M© for a wide range of the NS compact- 
ness if the BH has a prograde spin of a — 0.75. 
This amount of the disk mass for a high mass-ratio 
binary is hardly expected if the BH is nonspin- 
ning. This fact suggests that the formation of a 
BH-massive accretion disk system is a frequent out- 
come of the BH-NS binary merger with a prograde 
BH spin and may be encouraging for the merger 
scenario of a short-hard GRB. By contrast, the disk 
mass becomes very small if the BH has a retrograde 
spin. 

2. It is shown that some portion of the disrupted ma- 
terial can extend to > 400 km from the BH if the 
massive disk is formed. The maximum rest-mass 
density in the disk is larger for binaries with smaller 
values of Q because the ISCO radius is smaller for 
them. The extent of the disk could be large for 
a large value of Q. For such a remnant disk, the 
lifetime should be longer. 

3. The spin parameter of the remnant BH depends 
primarily on the spin parameter of the initial BH, 
a, and the mass ratio, Q. In particular, extrapola- 
tion of our results suggests that the merger of an 
extremely spinning BH and an irrotational NS does 
not form an overspinning BH. 

4. The gravitational waveform also depends strongly 
on the BH spin. The number of gravitational- wave 
cycles becomes larger for a prograde BH spin than 
that for a zero BH spin in the inspiral phase be- 
cause an additional repulsive force due to the spin- 
orbit interaction reduces gravitational-wave lumi- 
nosity and an approaching velocity of the binary. 
We found that the Taylor-T4 formula does not re- 
produce the phase evolution in the late inspiral 
phases accurately, especially when the mass ratio 
is large. 

5. In our previous work for nonspinning BH-NS bina- 
ries, the waveforms are classified simply into two 
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categories: when tidal disruption of the NS occurs, 
the waveform is composed of an inspiral waveform 
and a prompt shutdown at the tidal disruption. 
When tidal disruption does not occur, the wave- 
form is composed of inspiral and QNM waveforms. 
However, we find that the NS tidal disruption and 
the excitation of a QNM can occur simultaneously 
for binaries with a high mass ratio and a progradc 
BH spin. This is because the disrupted material 
cannot become axisymmetric before the prompt in- 
fall due to a larger BH areal radius for a larger value 
of Q. As a result, the material accretes onto the 
remnant BH coherently, and, therefore, the QNM 
of the remnant BH is excited, except for the case 
in which the extremely stiff EOS is adopted. 

6. The cutoff frequency of the gravitational- wave 
spectrum is correlated with the NS compactness 
in a clear manner when the NS is disrupted, and 
the BH spin modifies this correlation. The pro- 
grade BH spin decreases the cutoff frequency for 
fixed values of C and Q because the angular ve- 
locity at the tidal disruption becomes smaller than 
that for a — 0. The cutoff frequency is lower for a 
smaller value of C for fixed values of Q and a, as 
in the case of nonspinning BH-NS binaries, because 
the tidal effect is stronger and the disruption occurs 
at a more distant orbit. The BH spin also modifies 
the spectrum for the inspiral phase. Specifically, 
the spectrum amplitude for a given frequency in 
the inspiral phase becomes large when the BH has 
a prograde spin, and this is consistent with the PN 
estimation. Both the low cutoff frequency and large 
spectrum amplitude in the inspiral phase for a pro- 
grade BH spin are encouraging for gravitational- 
wave astronomy to become a tool to investigate 
the NS compactness and EOS. It is noteworthy 
that the BH-NS binary with a high mass ratio of 
Q > 5 is a more promising target for ground-based 
gravitational-wave detectors if the BH has a pro- 
grade spin and the NS tidal disruption occurs. 

Finally, we list several issues to be explored in the fu- 
ture. Piecewise polytropic EOSs with two pieces em- 
ployed in this paper are not accurate enough to model 
high-mass NSs with large central density of p max ^ 
10 15 g/cm 3 [59|. More detailed (piecewise polytropic) 
EOSs are necessary to calculate gravitational waves from 
a BH-rclatively massive NS binary merger for which the 
tidal deformation and disruption of the NS plays an im- 
portant role, i.e., BH-NS binaries with moderately large 
BH spins of a > 0.75. The implementation of detailed 
microphysics, such as a finite-temperature effect and a 
neutrino transport process, is essential even qualitatively 
to explore the evolution of the remnant BH-accretion disk 
system and to discuss the jet launch such that short- 
hard GRBs require. Recently, Sekiguchi developed a 
code to perform fully general relativistic simulation with 
the finite-temperature EOS and an approximate neu- 



trino emission (the so-called leakage scheme ) and succeed 
in simulating the stell ar c ore collapse 0, 1 1 OQj ] and the 
merger of binary NSs [101| . We plan to work on BH-NS 
binary mergers along these lines. 
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Appendix: Convergence of gravitational waves and 
the remnant disk mass 



This Appendix demonstrates that the convergence 
is approximately achieved for gravitational waves and 
masses of the remnant disks shown in Sec. IIVI Figure [25] 
shows the evolution of the orbital angular velocity deter- 
mined by Eq. Q32p and gravitational waveforms obtained 
with different grid resolutions for models HB-Q4M135a75 
and 2H-Q5M135a75. We perform an appropriate time 
shift in order to align the curves in the inspiral phase, and 
perform the rotation of + and x polarization modes of 
gravitational waveforms for N = 42 and 36. We note that 
the time to the merger, t mC rgcr, is systematically longer 
for finer grid resolutions because numerical dissipation 
of the angular momentum is smaller [2^, [3l|. In both 
cases, the evolution of the orbital angular velocity, f2(t), 
approximately agrees up to the merger, except for initial 
bursts associated with the junk radiation. The gravita- 
tional waveforms agree very well in the final ~ 4 orbits 
of the inspiral phase (t > 15 ms), the merger phase, and 
the ringdown phase if the QNM is excited. By contrast, 
gravitational waveforms in the initial ~ 2 orbits depend 
strongly on the grid resolutions, because the early inspi- 
ral phase is strongly affected by the dissipation of the 
junk radiation. We conclude that the convergence is ap- 
proximately achieved for gravitational waves in the late 
inspiral, merger, and ringdown phases. 

Figure [2U shows the evolution of the rest mass located 
outside the AH, M, 



r>TAM ' 



with different grid resolutions 
for models HB-Q4M135a75 and 2H-Q5M135a75. This 
figure shows that the convergence is also approximately 
achieved for the mass of the remnant disk. Quantita- 
tively, differences in M r>rAH at w 10 ms after the merger 
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arc w 2.9% and « 2.6% for models HB-Q4M135a75 and 
2H-Q5M135a75, respectively. If we assume the first-order 



convergence for M r>rAli , the errors in the values obtained 
for N = 50 runs are w 7.5% and ss 6.8% for models HB- 
Q4M135a75 and 2H-Q5M135a75, respectively. 
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FIG. 25. Comparisons of evolution of the orbital angular velocity and gravitational waveforms for models HB-Q4M135a75 
(top) and 2H-Q5M135a75 (bottom) with appropriate time shifts. 
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FIG. 26. Comparisons of evolution of the remnant disk masses 
for models HB-Q4M135a75 and 2H-Q5M135a75. 



magnitude of the pressure for the nuclear-density region. 
[63] B. D. Lackey, K. Kyutoku, M. Shibata, P. R. Brady, 

and J. L. Friedman, to be published. 
[64] M. Shibata, Prog. Theor. Phys. 96, 917 (1996). 



[65] P. Wiggins and D. Lai, Astrophys. J. 532, 530 (2000). 

[66] In this paper, "the ISCO radius" always represents "the 
ISCO radius in the Boyer-Lindquist coordinates," which 
is physical in the sense that it gives the proper cir- 
cumferential length for the equatorial circular orbit. It 
should be noted that the coordinate radius of the ISCO 
in our simulation is different from the Boyer-Lindquist 
one. 

[67] J. M. Bardeen, W. H. Press, and S. A. Teukolsky, As- 
trophys. J. 178, 347 (1972). 

[68] L. Blanchet, Living Rev. Relativity 9, 4 (2006). 

[69] T. Yamamoto, M. Shibata, and K. Taniguchi, Phys. 
Rev. D 78, 064054 (2008). 

[70] M. Shibata and T. Nakamura, Phys. Rev. D 52, 5428 
(1995). 

[71] T. W. Baumgarte and S. L. Shapiro, Phys. Rev. D 59, 
024007 (1998). 

[72] B. Briigmann, J. A. Gonzalez, M. Hannam, S. Husa, 
U. Sperhake, and W. Tichy, Phys. Rev. D 77, 024027 
(2008). 

[73] A. Kurganov and E. Tadmor, J. Comput. Phys. 160, 
241 (2000). 

[74] C. Reisswig and D. Pollney, arXiv: 1006. 1632. 

[75] M. Boyle, D. A. Brown, L. E. Kidder, A. H. Mroue, H. P. 

Pfeiffer, M. A. Scheel, G. B. Cook, and S. A. Teukolsky, 

Phys. Rev. D 76, 124038 (2007). 



33 



[76] L. Santamaria, F. Ohme, P. Ajith, B. Briigmann, 
N. Dorband, M. Hannam, S. Husa, P. Mosta, D. Pollney, 
C. Reisswig, E. L. Robinson, J. Seiler, and B. Krishnan, 
Phys. Rev. D 82, 064016 (2010). 

[77] L. E. Kidder, Phys. Rev. D 77, 044016 (2008). 

[78] T. Damour and A. Nagar, arXiv:0906.1769. 

[79] L. E. Kidder, Phys. Rev. D 52, 821 (1995). 

[80] L. G. Fishbone, Astrophys. J. 185, 43 (1973). 

[81] J.-A. Marck, Proc. Roy. Soc. London 385, 431 (1983). 

[82] M. Ishii, M. Shibata, and Y. Mino, Phys. Rev. D 71, 
044017 (2005). 

[83] J. E. McClintock and R. A. Remillard, in Compact 
Stellar X-ray Sources, edited by W. H. G. Lewin and 
M. van der Klis (Cambridge University Press, 2006). 

[84] K. Belczynski, T. Bulik, C. L. Fryer, A. Ruiter, 
F. Valsecchi, J. S. Vink, and J. R. Hurley, Astrophys. 
J. 714, 1217 (2010). 

[85] F. Foucart, M. D. Duez, L. E. Kidder, and S. A. Teukol- 
sky, Phys. Rev. D 83, 024005 (2011). 

[86] R. Penrose, in general relativity: an Einstein centenary 
survey, edited by S. W. Hawking and W. Israel (1979). 

[87] E. Berti, V. Cardoso, and A. O. Starinets, Class. Quan- 
tum Grav. 26, 163001 (2009). 

[88] K. Uryii, F. Limousin, J. L. Friedman, E. Gourgoulhon, 
and M. Shibata, Phys. Rev. Lett. 97, 171101 (2006). 

[89] K. Uryii, F. Limousin, J. L. Friedman, E. Gourgoulhon, 



and M. Shibata, Phys. Rev. D 80, 124004 (2009). 
[90] M. Saijo and T. Nakamura, Phys. Rev. Lett. 85, 2665 

(2000) . 

[91] M. Saijo and T. Nakamura, Phys. Rev. D 63, 064004 

(2001) . 

[92] M. Punturo et al, Class. Quantum Grav. 27, 194002 

(2010) . 

[93] L. Baiotti, T. Damour, B. Giacomazzo, A. Nagar, and 
L. Rezzolla, Phys. Rev. Lett. 105, 261101 (2010). 

[94] L. Baiotti, T. Damour, B. Giacomazzo, A. Nagar, and 
L. Rezzolla, Phys. Rev. D 84, 024017 (2011). 

[95] We refer to /dam as / cu t throughout in the previous work 
[3lj . In the present paper, we distinguish /dam from / cu t 
because the method for determining /d am is different 
from that for / cu t- 

[96] The relation between f cu t mp and C is different from the 
one obtained in Ref. [31[ due to the different definition 

Of /cut- 

[97] L. Rezzolla, Class. Quantum Grav. 26, 094023 (2009). 
[98] Y. Zlochower, M. Campanelli, and C. O. Lousto, Class. 

Quantum Grav. 28, 114015 (2011). 
[99] Y. Sekiguchi, Prog. Theor. Phys. 124, 331 (2010). 
[100] Y. Sekiguchi and M. Shibata, Astrophys. J. 737, 6 

(2011) . 

[101] Y. Sekiguchi, K. Kiuchi, K. Kyutoku, and M. Shibata, 
Phys. Rev. Lett. 107, 051102 (2011). 



